Dynamic Critical Behavior of the Chayes-Machta Algorithm 

for the Random-Cluster Model 



I. Two Dimensions 



Timothy M. Garoni 
Department of Mathematics and Statistics 
University of Melbourne 
Vic. 3010, AUSTRALIA 
T . GARONIOMS . UNIMELB . EDU . AU 

Giovanni Ossola 
Department of Physics 
New York City College of Technology 
300 Jay Street 
Brooklyn, NY 11201, USA 
GOSSOLASCITYTECH . CUNY . EDU 

Marco Polin 

Department of Applied Mathematics and Theoretical Physics 
University of Cambridge 
Wilberforce Road 
Cambridge CB3 OWA, UK 
M . POLINSDAMTP . CAM . AC . UK 

Alan D. Sokal* 
Department of Physics 
New York University 
4 Washington Place 
New York, NY 10003, USA 
SOKAL@NYU.EDU 

April 4, 2011 

Abstract 

We study, via Monte Carlo simulation, the dynamic critical behavior of the 
Chayes-Machta dynamics for the Fortuin-Kasteleyn random-cluster model, which 
generalizes the Swendsen-Wang dynamics for the g-state Potts ferromagnet to non- 
integer q > 1. We consider spatial dimension d = 2 and 1.25 < g < 4 in steps of 0.25, 
on lattices up to 1024^, and obtain estimates for the dynamic critical exponent ^cm- 
We present evidence that when 1 < q< 1.95 the Ossola-Sokal conjecture zcm > 
is violated, though we also present plausible fits compatible with this conjecture. 
We show that the Li-Sokal bound zqm > a/u \s close to being sharp over the entire 
range 1 < (7 < 4, but is probably non-sharp by a power. As a byproduct of our work, 
we also obtain evidence concerning the corrections to scaling in static observables. 

Key Words: Random-cluster model. Potts model, Chayes-Machta algorithm, Swendsen- 
Wang algorithm, cluster algorithm, dynamic critical behavior. 
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1 Introduction 



Since nontrivial models in statistical mechanics are rarely exactly solvable, Monte 
Carlo (MC) simulations have become a standard tool for obtaining information on phase 
diagrams and critical exponents [lH3]. Unfortunately, MC simulations typically suffer 
from severe critical slowing- down |11[5], so that the computational efficiency tends rapidly 
to zero as the critical point is approached. More precisely, the autocorrelation (relaxation) 
time r diverges in the critical limit, most often like r ~ where ^ is the spatial correlation 
length. The dynamic critical exponent z depends on both the model being investigated 
and the MC algorithm being used. For local algorithms one typically has z ^ 2. 

An important advance was made in 1987 with the invention of the Swendsen-Wang 
(SW) cluster algorithm [B] for simulating the ferromagnetic g-state Potts model [THH] at 
integer q > 2. The SW algorithm is based on passing back and forth between the Potts 
spin representation and the Fortuin-Kasteleyn (FK) bond representation [TUHT^ . More 
precisely, one introduces a joint probability distribution [H] of spin and bond variables, 
whose marginal on the spins (integrating out the bonds) is the Potts spin model and 
whose marginal on the bonds (integrating out the spins) is the FK random-cluster model 
[U]; one then updates this joint distribution by alternately applying the two conditional 
distributions (see Section fLTl below for details). 

Since a local move in one set of variables can have highly nonlocal effects in the 
other, it is not surprising that the SW algorithm might have less critical slowing-down 
than the conventional local algorithms. And this is in fact the case: although the SW 
algorithm does not eliminate critical slowing-down, it does radically reduce it compared 
to local algorithms. Much effort has therefore been devoted, for both theoretical and 
practical reasons, to understanding the dynamic critical behavior of the SW algorithm as 
a function of the spatial dimension d and the number q of Potts spin states. The best 
information on the dynamic critical exponent 2;sw prior to the present work is summarized 
in Table [TJ Unfortunately, it is very difficult to develop a physical understanding from the 
small number of non-trivial "data points" at our disposal: second-order non-mean-field 
transitions occur only for (rf, q) = (2, 2), (2, 3), (2, 4), (3, 2) and (4, 2)0 

A further advance was made in 1998 by Chayes and Machta (CM) [25], who devised 
a cluster algorithm for simulating the FK random-cluster model at any real value q > 1. 
The idea behind the CM algorithm is very similar to that of SW, but now one starts 
with the random-cluster (bond) measure and introduces auxiliary color (spin) variables 
on the sites, so that we again have a joint model of spins and bonds, analogous to the 
one employed in the SW algorithm. Although the marginal measure on the spins no 
longer has any obvious physical interpretation when q is noninteger, the joint measure 
can still be used to construct an efficient cluster algorithm, by alternately applying the 
conditional distributions exactly as in standard SW. The CM algorithm thus generalizes 
the SW algorithm and in fact reduces to (a slight variant of) it when q is an integer: 
see Section 12.21 for details. Indeed, the CM algorithm can be thought of as a "natural" 



^ A continuous (second-order) transition occurs also in the Ising {q = 2) model in dimensions d > 4, 
but here the static behavior is mean-field. One expects the dynamic critical exponents likewise to be 
dimension-independent for d > 4 (with possible multiplicative logarithmic corrections at o? = 4) . 
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Table 1: Best estimates of the dynamic critical exponent z for the Swendsen-Wang 
algorithm prior to the present work. Estimates are taken from [16] for d = 2, q = 2; [T7] 
for = 2, g = 3; for d = 2, g = 4; [20] for d = 3, g = 2; and ^ot d = A, 

q = 2. Error bars are one standard deviation, and include statistical errors only. 



interpolation of the SW algorithm to noninteger q (though unfortunately only for q > 1). 
Thus, by using the CM algorithm we can study the dynamic critical behavior of the SW- 
CM dynamic universality class as a function of the continuous variable q throughout the 
range 1 < q < qdC)-, where gc(£) is the maximum q for which the transition is second- 
order on the given lattice Cn This vastly enhances our ability to make theoretical sense 
of the numerical results. 

In the present paper we perform high-precision MC simulations of the random-cluster 
model on the square lattice {d = 2), for 1.25 < g < 4 in steps of 0.25, using the CM 
algorithm. Our main goal is to gain better insight into the SW-CM dynamic universality 
class in two dimensions by studying the behavior as a function of the continuous variable q. 
We estimate numerically several dynamic critical exponents for each value of q, and we 
attempt to understand their behavior as a function of q. As a byproduct we also obtain 
new information on corrections to scaling in the static quantities. In the companion paper 
[27] we carry out an analogous study for the simple-cubic lattice (rf = 3) at g = 1.5, 1.8, 2.2; 
and in a forthcoming paper [28] we will analyze the case of the complete graph (Curie- 
Weiss model) jfl 

One advantage of considering c? = 2 is of course is that we have available a remarkable 
amount of information concerning the static behavior of the random-cluster model as a 
function of g. This information includes the exact location [30l|31] of the transition point 
on the square lattice, namely pc = ^/q/i)- + y/q), and the knowledge that the transition 
is second-order for < g < 4 and first-order for g > 4. Furthermore, the leading and 
next-to-leading thermal and magnetic critical exponents are known exactly (though non- 

^ We stress that qdC) is not necessarily the same for all lattices of a given dimension d; the first-order 
or second-order nature of the transition is a non-universal question. See [26j for further discussion. For 
the standard two-dimensional lattices (square, triangular and hexagonal) we have qc{C) — 4. 

A brief summary of the preliminary results in this trio of papers has appeared in letter form [29] . 
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= 2yHi-d = (1.10 

(Here Ai is the leading correction-to-scaling exponent, defined e.g. by corrections ~ L^^^ 
in finite volume at criticality; and dp is the cluster fractal dimension.) The numerical 
values of these critical exponents for the values of q employed in our simulations are 
collected for reference in Table [2l 

A key theoretical result concerning the SW dynamics is the Li-Sokal bound [35j, which 
states that 

T-int,Ar, T-exp,Ar > COUSt X Cr (l-H) 

where Tintx ^^"^ ^e^v,M are, respectively, the integrated and exponential autocorrelation 
times for the observable M = number of occupied bonds, and Ch is the specific heat. (For 
the precise definitions of Tint and Tgxp, see Section ISTTl ) It follows from this lower bound 



See |191 Appendix A.l] for further discussion and references. Note that there is a typographical error 
inequation (A. 10) of [19] . which should read ^ = {[2{s — r) + sx]'^ — x^)/{8[{2 + x)]. We remark that the 
same formulae for 4 < g < 6 give the exponents yri, 2/T2, y//i, yH2 of the tricritical Potts model |32 | I33 ) . 
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that the SW algorithm cannot completely eliminate critical slowing-down if the specific 
heat is divergent at criticality. In particular, (11. lip implies the lower bound 



Zint,M^ ^exp > a/l^ (1-12) 

for the dynamic critical exponents of the SW algorithm, where a/u is the static exponent 
for the specific heat. Now, simple "Fortuin-Kasteleyn identities" show that var(A/')/K is 
a specific-heat-like quantity and provides a natural continuation of the notion of specific 
heat to random-cluster models. It turns out that the Li-Sokal bound (11.111) with this 
definition of specific heat can be easily extended from the SW to the CM algorithm at 
arbitrary real g > 1, as we show in j271 Appendix A]. 

The physical mechanism underlying the Li-Sokal proof is the slow evolution of A/", 
which is an "energy-like" observable. Previous studies in two dimensions for q = 2,3,4 
have shown empirically that the Li-Sokal bound is very close to being sharp for all three 
values of q [HI Section 6] [T7HT9]; and we will show here that this is the case also for 
noninteger q throughout the range 1 < g < 4 (see Sections 17.31 and I7.4p . For the three- 
dimensional Ising model, by contrast, the Li-Sokal bound is very far from sharp |2Uj . 
Therefore there must be another mechanism, beyond the one captured in the Li-Sokal 
proof, that is principally responsible for the critical slowing-down in three dimensions. A 
few years ago, Ossola and Sokal |20] suggested that this as-yet-not-understood mechanism 
causing slowness might perhaps be somehow related to the typical size of the largest 
cluster. More specifically, they conjectured that perhaps the SW dynamics for any Potts 
ferromagnet in any dimension satisfies 

Tsw > const X — (1-13) 

where Ci is the expected number of sites in the largest cluster. If true, this would imply 
the critical-exponent inequality 

^sw > (1.14) 

where P/u is the static exponent for the magnetization. Indeed, Coddington and Baillie 
[23] had earlier suggested that the inequality (I1.14p holds as an equality for the Ising 
models in dimensions d = 2,3,4. The numerical results reported by Ossola and Sokal 
[20] for the three-dimensional Ising model are consistent with this conjectured equality, 
though they also present plausible fits consistent with < P/i^ (i-e. violation of their 
conjectured inequality). Finally, an analytical treatment [23] of the SW dynamics for the 
Ising model on the complete graph suggests that 2;sw = P/^ in this case also (namely, 
^sw = P/t^ = 1)- 

Alas, by considering the Ossola-Sokal conjecture within the more general framework 
of the Chayes-Machta dynamics, one can see at a glance that it is probably false! To start 
with, for g = 1 on any lattice, the Swendsen-Wang-Chayes-Machta algorithm reduces to 
independent sampling for independent bond percolation, so that zqm = 0; but in general 
/3/z/ > for percolation, which shows that the Ossola-Sokal conjecture (I1.14p is false for 
q = I. Moreover, it is reasonable to believe that zcm is a continuous function of q: if so, 
then for q slightly greater than 1 one has either zqm identically zero or else ^cm very close 
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to zero; and in either case one would have zqm < /9/za so that the Ossola-Sokal conjecture 
would be violated also for q in some interval above l|j Our numerical results in this paper 
(see Section [775]) suggest that in d = 2 the Ossola-Sokal conjecture fails for 1 < g < 1.95; 
indeed, for 1 < g < 1.6 we apparently have zqm = exactly, while of course /3/z/ > 0. (But 
see Section [775] for an alternative fit that is compatible with the Ossola-Sokal conjecture.) 
Likewise, the numerical results for = 3 to be presented in [27] for q = 1.5 and q = 1.8 
strongly suggest that the Ossola-Sokal conjecture is violated for these values of q and 
hence presumably for all g < 2, but that it holds as a strict inequality zqm > P/t^ when 
q = 2.2 and hence presumably for all q in the range 2 < g < gc(Z^). Indeed, in ci = 3 it 
is found that zqm is a strongly increasing function of q, while P/u varies very slowly with 
q; the two curves presumably cross somewhere near q = 2. Finally, in our forthcoming 
paper [28] we will extend the method of [2l] to cover the CM dynamics on the complete 
graph for all q in the range 1 < g < 2. We find that there are two markedly different 
behaviors depending on g: for 1 < g < 2 we have zqm = < fS/u = 1, while for g = 2 we 
have zcM = (3/^ = 1. In particular, the Ossola-Sokal conjecture fll.Mp is violated when 

1 < g < 2 for the CM dynamics on the complete graph. 

It therefore seems that the mechanism causing slowness of SW-CM in three dimensions 
is not related in any simple way to the size of the largest cluster, contrary to the intuition 
of Ossola and Sokalj^ And so we are back to square one: no one seems to have the slightest 
idea what is the physical mechanism dominating the critical slowing-down of the SW-CM 
dynamics in dimensions d > 3. 

The present paper is organized as follows: In Section [2] we briefly review the Swendsen- 
Wang algorithm and then present a simple explanation of the Chayes-Machta algorithm. 
In Section [3] we define the observables and autocorrelation times that we measured in 
our simulations, while in Section H] we discuss our methods of statistical data analysis. In 
Section[5]we summarize the characteristics of our MC simulations. In Section[6]we analyze 
our numerical data for the static observables, with an emphasis on detecting corrections 
to scaling and estimating their exponents and amplitudes. In Section [7] we analyze our 
numerical data for the dynamic quantities (i.e. autocorrelation times), with an emphasis 
on estimating the dynamic critical exponent ^cm as a function of g; in particular we test 
the Ossola-Sokal conjecture and the sharpness of the Li-Sokal bound. Finally, in Section[8] 
we briefly discuss our results and the prospects for future work. 

2 The Swendsen— Wang and Chayes-Machta algorithms 

In this section we describe a family of algorithms for simulating the FK random-cluster 
model at any real g > 1, which generalize slightly the original algorithms introduced 
by Chayes and Machta [25] . We begin (Section 12. ip by reviewing the Swendsen- Wang 



^ It is of course conceivable that zcm is discontinuous at q ~ 1: for instance, one might have zcm = 13 /i^ 
exactly for all q near 1, but with an amplitude that vanishes as g J, 1. In Section [7?5] we will test this 
scenario against our data. 

^ Or at least, no such relation seems to hold for all q. It is still conceivable that such a relation might 
hold for q = 2 only, or for all q > 2. 
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[6] algorithm for simulating the ferromagnetic g-state Potts model at integer q > 2, 
with emphasis on the perspective afforded by the joint probability distribution of spins 
and bonds (sometimes known as the Edwards-Sokal coupling [H]). We then present 
(Section I2.2p our version of the Chayes-Machta algorithm and a proof of its validity. For 
an alternative (and rather more general) presentation of the Chayes-Machta algorithm, 
see [27]. 

2.1 Swendsen— Wang algorithm 

Let G = (V, E) be a finite graph with vertex set V and edge set E, and let q be an 
integer > 2. The g-state Potts model on G with nearest-neighbor couplings f3 = {f3e}eeE 
is defined by the Gibbs measure 

Hfi^q{cr) (X exp I ^/3e5e(cr) I (2.1a) 

eeE 

for cr G {1, 2, . . . , q}^ , where Pe = 1 — e"^"" for each e & E and 

W = S^.,<ry for e = (xy) . (2.2) 

One rather natural way to elucidate the Swendsen-Wang (SW) algorithm is via the 
Edwards-Sokal coupling [13] of the Potts and random-cluster models. Recall first that the 
random-cluster model [15], introduced by Fortuin and Kasteleyn [TOl - fTS] . is a correlated 
bond-percolation model defined on a finite graph G = {V, E) for parameter g > (not 
necessarily an integer) and edge probabilities p = {pe}e€E, by the probability measure 

0p,,(n) (X g'=(")5p(n), (2.3) 

for n G {0, l}'^, where 

BA^) =11^^11 (1-^^) (2.4) 

e£E: ne=l e£E: 71^=0 

and k{n) is the number of connected components ("clusters") in the spanning subgraph 
whose edges are those with ng = 1 (we call these "components of n" for short). The Potts 
and random-cluster models are intimately related, in the sense that correlation functions 
of the Potts spins can be expressed in terms of connectivity functions of the random- 
cluster model. The Edwards-Sokal coupling is a joint measure on the space of spin and 
bond configurations, defined by 

Hp,q{n,(T) oc Yl[il-Pe)5n„0 + Pe5n,,l5e{a)] . (2.5) 
e£E 
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It is easy to show j5|[T^ that the marginal measure '^nt^p,g('^' ') spins is the Potts 

measure fl2.ip . and the marginal measure J^a- /^p,g('' ''') bonds is the random-cluster 

measure fl2.3l) . Explicit forms for the conditional measures are also easily computed: 



where 



/ip,,(o-|n) = A(n,o-)g-'=(-) (2.6) 
/ip,,(n|cr) = A(n,cr)5p(n) (2.7) 



, ^ if 5Jcr) = 1 for all e with n^, = 1 , . 

^l^''^) = <!. J . (2-8) 
otherwise 



In words, A(n, cr) is the indicator for the event that "we draw occupied bonds only 
between vertices with the same spin value". 

The SW algorithm simulates the Edwards-Sokal measure (12. 5p . and hence both the 
Potts and random-cluster models, by alternately updating the spins conditioned on the 
bonds using (12. 6p and the bonds conditioned on the spins using (12. 7p . In words, the 
spin- updating rule (12. 6p states that, independently for each cluster, we choose a new spin 
value uniformly at random from the q possible choices and assign this spin value to all 
vertices in the cluster; while the bond-updating rule (12.70 states that on each subgraph 
induced by the vertices of a given spin value we update the bonds via independent bond 
percolation. 



2.2 Chayes— Machta algorithm 

Consider the random-cluster model defined on a finite graph G = {V,E) by (12. 3p . 
Recall that to arrive at the SW algorithm for the Potts model, one introduces auxiliary 
bond variables and then considers the joint model of Potts spins on the vertices and 
auxiliary bond variables on the edges. The SW algorithm updates this joint model by 
alternately applying its two conditional measures (bond variables given spins and spin 
variables given bonds). Here we do the reverse: starting from the random-cluster model 
(12. 3p . we introduce auxiliary color variables on the vertices of G and then consider the 
joint model of our original bond variables and the new auxiliary color variables; we will 
then update this joint model by alternately applying its conditional measures. To this 
end, let us introduce color variables on the sites taking values in some finite set S, i.e. 
cr e . For each a E S, we choose a number qa > 0, and write q = {qa}aes- We can 
now introduce the coupled probability measure 

/xp,q(n, ct) oc A(n, cr) B^{n) J] ^''^^"'"^ • (2-9) 

Here A(n, cr) is again given by (12.81) but now with with cr G , i.e. it is the indicator 
of the event "we draw occupied bonds only between vertices of the same color" ; and we 
denote by ka{cr,n) the number of components of n that are colored a by cr. Note that 
ka{cT,n) is well-defined whenever A(n,cr) ^ 0. From the definition (12.90 we can observe 
the following: 
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1) The marginal measure Yla- f^p,<ii' ' ^) °^ bond variables is the random-cluster 
measure with parameter q = 



2) The conditional measure of cr given n is as follows: Independently for each compo- 
nent of n, randomly choose a value a & S with probability Qa/Q and impose it on 
all vertices in that component. 

3) Every coloring cr defines a partition of the vertex set into V = {J^es where Va 
is the set of all sites v E V colored a by cr. The conditional measure of n given cr 
is then simply the product over a G 5* of random-cluster measures with parameter 
Qa on the induced subgraphs All bonds not lying within a single G[Va] are 
forced to be unoccupied. 

Now suppose that for each a E S and each induced subgraph H C G, we have a 
valid Monte Carlo algorithm for updating the random-cluster model with parameter qa 
on H. We can then apply the following algorithm to simulate the joint model f l2.9p . and 
therefore, via Observation 1, the random-cluster model with parameter q = Xlags^'a- 

Algorithm 1 (Random-cluster algorithm). 

1. Given a bond configuration, choose a new color configuration using the conditional 
measure described in Observation 2. 

2. Given the new color configuration, for each a E S update the random- cluster model 
on G\ya\ with parameter q^, using the given algorithm. 

Note that for any a E S, one valid Monte Carlo update is to do nothing (i.e. perform 
the identity operation), since detailed balance is satisfied trivially. Furthermore, if q^ 
happens to equal 1, one can use a Bernoulli update, i.e. erase the bonds on G\Va\ and 
choose new ones via independent bond percolation. But if one has some other way of 
updating a random-cluster measure with parameter then that is fine too. A number 
of special cases of Algorithm [1] are now clear: 

1. If q is an integer, S" = {1, . . . , g}, ga = 1 for all a, and one uses Bernoulli updates 
for every a, then Algorithm [1] is simply the standard SW algorithm. 

2. If q is an integer, S" = {1, . . . , g}, = 1 for all a, and one uses Bernoulli updates 
for 1 < a < k and do-nothing updates for k + 1 < a < q (where k is any integer 
satisfying 1 < A; < g), then Algorithm [1] is a variant of the SW algorithm in which 
only the first k colors are "active" . 

3. If g > 1, 5* = {active, inactive}, qactive = 1, QinacUve = Q — ^, and one uses a Bernoulli 
update for a = active and a do-nothing update for a = inactive, then Algorithm [1] 
is the original Chayes-Machta algorithm presented in [25j. 

4. If g > A; > 1 with k an integer, 5* = {1, . . . ,k, inactive}, g^ = 1 for 1 < a < /c, 
^inactive = 1 — k, and one uses a Bernoulli update for 1 < a < A; and a do-nothing 
update for a = inactive, this is the generalized algorithm presented in [251 PP- 482- 
483]. 
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Note that since each of these cases involves performing Bernoulh updates, they all 
require q > 1. For a given value of g > 1, Case 4 provides a family of algorithms, indexed 
by the integer k, for simulating the random-cluster model. We call this "the Chayes- 
Machta algorithm with k active colors" . Any choice of 1 < k < \_q\ is legitimate, though 
it seems reasonable to expect that k = \_q\ should be the most efficient. Indeed, our 
numerical simulations suggest that in practice the autocorrelation time is approximately 
proportional to 1/k (see Section [711 below) . 



3 Observables and autocorrelation times 

In this section we recall the definitions of the various autocorrelation times (Sec- 
tion [3?T]) and then list the observables that we measured in our simulations (Section l3.2p . 



3.1 Autocorrelation functions and autocorrelation times 

Consider an observable O in the random-cluster model, i.e. a real-valued function of 
n G {0, l}'^. Then a realization of the Chayes-Machta (CM) Markov chain gives rise to 
a time series 0{t), where each unit of time corresponds to one step of the CM algorithm. 
The autocovariance function of O is defined to be 

Cooit) = (O(O)O(t)) - (3.1) 

where the expectation is taken in equilibrium. The normalized autocorrelation function 
of O is then 

Poo(t) = (3.2) 
From pooit) we define the integrated autocorrelation time as 

Tint,o = 2 Poo{t) (3.3) 

t=~oo 

and the exponential autocorrelation time as 

Texp.o = limsup- ^^Wtt- (3-4) 

\t\->oo log pAfAf{t) 

Finally, the exponential autocorrelation time of the system is defined as 

Texp = SUpTexp.O, (3.5) 
O 

where the supremum is taken over all observables O. This autocorrelation time thus 
measures the decay rate of the slowest mode of the system. All observables that are not 
orthogonal to this slowest mode satisfy Tg^p^o = ^exp- 

It is important to remember that there is not just one autocorrelation time, but many: 
namely Texp as well as Tint,© for each O. In all but the most trivial Markov chains these 
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autocorrelation times are not equal. Correspondingly, there are many dynamic critical 
exponents: namely z^^p as well as ^mt.o for each O. These exponents may in some cases 
be equal, but they need not be; this is a detailed dynamical question, and the answer will 
vary from algorithm to algorithm and model to model. 

More information on the principles of Markov-chain Monte Carlo and the relations 
between these autocorrelation times can be found in [5J. 



3.2 Observables to be measured 

We now take the graph G to be the dimensional simple-hypercubic lattice of linear 
size L with periodic boundary conditions; we denote this graph by in order to em- 
phasize that addition of coordinates is always taken mod L. When the precise form of 
the graph G is unimportant to our definitions, we denote the vertex set and edge set as 
simply V{G) and E{G), respectively. 

For xi, . . . ,Xm G V{G), not necessarily all distinct, we denote by 'yxi---xm the indicator 
for the event that the vertices xi, . . . , Xm are all in the same cluster, i.e. 

, , 1 1 if Xi ^ a;2 ^ ^ ■ ■ ■ ^ in configuration n 
lx,-xjn) = I ' _ ' (3.6) 

I otlierwise 

and X y denotes that x is connected to y by at least one path of occupied bonds. In 
particular, T{x,y) = {'jxy) is the "two-point connectivity function", i.e. the probability 
that X is connected to y. On we write t{x) = t{x,0), and translational invariance 
gives us t{x, y) = t{x — y) for all x,y E Z^. 

We measured the following observables in our MC simulations: 



The number of occupied bonds 



Af = J2 rie (3.7) 

eeEiG) 



The nearest-neighbor connectivity (which is an energy-like observable [T7] ) 

{xy)&E{G) 



The cluster-size moments 

Sm = 1^1™ (3.9a) 

C&K{n) 

Yl (3.9b) 

Xi,...,XmeV{G) 

where the size \C\ of a cluster means the number of vertices. In this work we 
measured Sm for m = 2,4, 6, 8; 
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The size Ci of the ith largest cluster. In this work we measured Ci for i = 1,2, 3; 

An observable used to compute the Fourier transform of T{x,y) evaluated at the 
smallest non-zero momentum {2it/L, 0, . . . , 0): 



1 



7xy e 



2TTi{xj-yj)/L 



d 



x&C 



(3.10a) 
(3.10b) 



i=l C&K{n) 

From these observables we computed the following quantities: 
• The bond density 



N 



(AO 
B 



(3.11) 



where B = \E\ is the total number of edges in the graph G (i.e. B = dL'^ for 

G = ziy, 

The connectivity density 



E' 



B 



(3.12) 



The specific heat. There are a number of sensible definitions for specific heat, and 
we considered the following (see Remark 13.11 below) : 



cS'^ = 52ic'H'-i^-p)m 



-var(Ar) 



d 



(1) 



P 



C 



(3) 
H 



C 



(2) 



(3.13) 
(3.14) 

(3.15) 



(g-l)2 ^ 

(2) 

In this paper we used Ch = Ch in agreement with 
The expected size of the ith largest cluster, 

a = (3.16) 
We expect that in the critical region Ci ~ L'^-P/^ as L — )■ oo, for any fixed i; 
The susceptibility (see Remark [3^2] below) 

{S2) 



X 



(3.17) 



In the random-cluster model, x is the mean size of the cluster containing any spec- 
ified point. 
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The Fourier transform of the two-point connectivity function t{x), evaluated at the 
smallest nonzero momentum, {27t/L,0, ...0): 

F = r{p) = ^3.18) 

p={27r/L,0,...0) 

The finite-size second- moment correlation length (see Remark 13.31 below): 

1) (3.19) 



2sin(7r/L) \F 



Remark 3.1. The definition C^^ of the specific heat seems quite natural from the perspec- 
tive of the random-cluster model, in which we view the bond variables as fundamental; 
in particular, it arises in a very direct way in the proof of the Li-Sokal bound for the 
CM algorithm [271 Appendix]. The definitions C^^ and by contrast, are designed to 
reduce to more familiar expressions in terms of the energy of the Potts spin model when 
q is an integer. Specifically, if q is an integer we have 

Cj?) = -varl >J 5.^,.^ 1 (3.20) 




CP = -varl a, -ay] (3.21) 

\{xy)GE / 

where ax G {1, 2, . . . , g} and cr^. G M^~^ is a Potts spin in the hypertetrahedral represen- 
tation, so that 

<^x-<^v = ^ """^'""^ . (3.22) 

q — I 

Remark 3.2. We note that for q G {2, 3, . . .}, the Fortuin-Kasteleyn identity (cr^,. ■ cTy) = 
{jxy) shows that {M"^) = {S2), where A^^ is the squared magnetization, and hence x is 
simply equal to the magnetic susceptibility. 

Remark 3.3. In terms of r(p) the definition (13.1 9p has the general form 

1 ( ^(0>--->0) (3 23) 

^ 2sin(7r/L) Vr(27r/L,0,...,0) J ^ ' 

and applies equally well to the integer q case where we interpret r(x) as the spin-spin 
correlation function. A very nice account of the origin of the definition fl3.23p of the 
finite-size correlation length, relating it to the thermodynamic second-moment correlation 
length is given in [36] Part III]; see also [371138] . 

The prime on J-"' is employed to distinguish it from another estimator for F, often 
denoted J-", which can be defined in terms of Potts spins when g G {2, 3, . . .}: see e.g. [T7] . 
We note that J-'' could profitably be used even in Swendsen-Wang simulations of the Potts 
model, as an alternative to J-". 
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Remark 3.4. For any g > and any (xy) G E{G) we have the identity 

q{n.y) = p[(g-l)(7.,) + l]. (3.24) 

This is most easily seen by first considering the Edwards-Sokal joint measure when q is 
an integer > 2, and then combining the identities (uxy) = p{Sa^,ay) and {'jxy) = {(^x ■ CTy) 
together with f l3.22p . See e.g. [S]. This clearly proves f l3.24p for all q G {2,3, . . .}; and 
since both sides are rational functions of g, it follows that the equality must hold for all 
g G C. Summing over all edges we then obtain 

{U) = p^^{£') + -B. (3.25) 
q q 

As a consistency check on the correctness of our simulations we numerically tested the 
identity (I3.25P to high precision by measuring the observable 

U - p'^^E' - ^B . (3.26) 
q q 

Clearly this observable should have mean zero, and in all our simulations this was observed 
to be the case within statistical errors. 



4 Statistical analysis of the Monte Carlo data 

Consider a generic observable O with expectation (C) = jio- Suppose that, after 
equilibration, we have a sequence of T Monte Carlo measurements of C, denoted {Cf}^^. 
Since the Monte Carlo process is assumed equilibrated, the {Ot\^^i are a sample from 
a stationary stochastic process (also called a stationary time series). In this section we 
recall the basic principles of statistical time-series analysis, and describe the standard 
estimators for various quantities associated with O. For more details, see e.g. [5|[39]. 

The natural estimator for the expectation /zq is the sample mean 

= ^Y.^t. (4.1) 

t=i 

This is an unbiased estimator and has a variance 

T 



var(0) = l-^Y.Coo{t-s) (4.2a) 

S,t=l 

t=-(T-l) ^ ^ 

~ ^2ri„t,oCoo(0) for T > ri„t,o • (4.2c) 

This implies that the variance of (9 is a factor 2 Ti^t,o larger than it would be if the 
measurements were uncorrelated. Therefore, in order to obtain accurate error bars for 
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the estimator of the static quantity j^o, we must obtain an accurate estimate of the 
dynamic quantity Ti^t,o- 

The natural estimator for the autocovariance function is 

. T-\t\ 

Cooit) = 5^(a-/io)(a+*-/io) (4.3) 

if the expectation no is known, and 



1 _ _ 

Cooit) = j—-^ - ^)(^^+* - ^) (4.4) 



if /io is unknown. The estimator Coo{i) is unbiased, and the bias of Cooi^) is of order 

1/T. To leading order for T ^ Tint,© the covariance matrices of Coo and Coo are equal 
and it can be shown [39lll0l that 



cov{Cooit), Cooiu)) = Tj; Y [Cooi-s] Cooi-s + u - t) + Cooi-s + u) Cooi-s - 

+ K{t,s,s + u)] + o{l/T) , 
where t, m > and k is the connected 4-point autocorrelation function 
K{r,S,t) = {{Oi - Ho)iOi+r - fJ'0){0i+s - fio)iOi+t - fJ'o)) 

- Coo{r)Coo(t - s) - Coo{s)Coo(t - r) - Coo(t)Coo{s - r). 
Similarly, the natural estimator for the autocorrelation function is 



(4.6) 



- Coo{t) 



if lio is known, and 



= (4.8) 
Coo(O) 

if /io is unknown. Both the estimators 'poo{t) and ^ooi^) have bias of order 1/T. The 
covariance matrices of ^oo and 'poo are the same to leading order for large T. If the 
process is Gaussian, this covariance matrix is given in the large T limit by [39j 

^ oo 

co\{pooit)^Poo{u)) = 7^ Y [poo{.s)poo{.s + t-u)+ poois + u)poois - t) 

+ 2 pooit)pooiu)pl,ais) - 2 poo{t)poo{s)poo{s - u) 
-2poo{u)poo{s)poo{s -t)] + o(l/r) . 

If the process is not Gaussian, then there are additional terms proportional to the fourth 
cumulant K,{s,t,t — u). The simplest assumption is to assume the stochastic process to 
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be "not too far from Gaussian" and simply drop all the terms involving k. This is what 
we will do. If this assumption is not justified, then we are introducing a bias into the 
estimate of cov (pooit), Poo (u)). 

Finally, we take the estimator for the integrated autocorrelation time to be [H] 

1 ^ 

t=-M 

if /xo is known, or the analogous object defined in terms of is fio is unknown; here 
M (< T — 1) is a suitably chosen number, which we call the window width. One's first 
thought might be to use all the data, i.e. take M = T — 1; but this turns out to be a 
very bad estimator, which has a variance of order 1 even as T — > oo. Loosely speaking 
(see e.g. [39lllT]), this is because the terms in pooit) with large t (namely, t ^ Tcxp,©) 
have variance of order 1/T that does not vanish as t grows, cf. (14. 9p . and there are of 
order T of them. These terms thus contribute much "noise" but very little "signal" , since 
Pooit) is tiny when t ^ Texp,©- To obtain a good estimator, we should instead choose the 
window width M to be large enough so that we do not lose too much "signal", i.e. pooif) 
is tiny for all t > M, but not too much larger than this. In general, a window of width 
M creates a bias given by 

bias(fint,o) = -^Yl P^^^^^ + ^(V^) • (4.11) 

\t\>M 

The variance of the estimator Tij^t,o can be computed from (14. 9p : assuming that Ti^t,o ^ 
M -C T, the final result is 01] 

, 2(2M + 1) . 
var(rint,o) ^ '-tL,o ■ (4-12) 

To choose the window width M we used the automatic windowing algorithm intro- 
duced in [H], in which one sets 

M = min{m G Z : m > CTint,o{fn)} (4-13) 

where c is a suitably chosen constant. If the normalized autocorrelation function is ap- 
proximately a pure exponential, then a choice in the range c ~ 6 — 8 is reasonable. 
Indeed [17], if we take poo{i) = e"'*'/'^ and minimize the mean-square error 

MSE(fint,o) = bias(fint,o)^ + var(fint,o) (4.14) 
using (14. lip / (I4.12p . we find that the optimal window width is 

Mopt = 2^°s(^) -1- (4.15) 

For n/r ^ 10^ 10^ 10^ 10^ 10^ with r > 1, we have M^pt/r ^ 4.26, 5.41, 6.56, 7.71, 8.86, 
respectively. 
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In this paper we chose c = 6 for the observable £' , which has the slowest-decaying 
autocorrelation function of all the observables we measured and whose autocorrelation 
function is in fact very close to a pure exponential. We then used the same window 
width M, computed by applying the automatic windowing algorithm to £' , for all the 
other observables. We think that this procedure is preferable to applying the automatic 
windowing algorithm directly to the other observables, since some of the latter have an 
autocorrelation decay that is far from a pure exponential (this is especially so for J-"', C2 
and C3). 

Remark 4.1. If . . . , C^'")} is a family of observables, we estimate composite observ- 

ables of the form /((C^^^), . . . , (C^™))) by /(O^, . . . , OM) and approximate the corre- 
sponding variance by assuming the fluctuations from the mean are small, so that 



where 



var/(CW,...,CM) ^ ^^varZ, (4.16) 



5^a.O? (4.17) 



i=l 



^ ^^i{0,),...,{Om)) (4.18) 



and g = g{{Oi) , . . . , (Om)) is just a convenient common factor of the derivatives of /, 
which may be simply 1 in practice. Obviously in practice we compute using Oi in place 
of {Oi). We used this method to compute the estimate and standard error of ^. 

In principle we can also estimate the specific heat using this procedure, but in practice 
it turns out that a better estimator for the specific heat is provided by 

Ch = ^i^-m'. (4.19) 

The mean of KW} is clearly CjJ^ [cf. (I37[3|) ]. and its variance determines the error bar 
on the resulting estimate of C^'*. We prefer this estimator since it avoids computing the 
quantity A/"^ — A/° which can be a numerically perilous object when jV^ and AT are close 
in value. 



Remark 4.2. For the larger values of L we actually performed a number of independent 
runs rather than one long run. The best estimate of each autocorrelation function poo{t) 
was then computed, for each fixed t, by averaging the ^ooi^) from the individual runs, 
with weights proportional to the run lengths. The windowing procedure described above 
was then applied to this best estimate of pooii) in order to obtain the final value for the 
integrated autocorrelation time. 



18 



5 Description of the simulations 



We implemented the Chayes-Machta (CM) algorithm for the random-cluster model on 
an L X L square lattice with periodic boundary conditions. We performed all our runs at 
the exact critical temperature pc = y/q/ (1 + y/q). We studied lattice sizes 16 < L < 1024 
in powers of 2, and parameters 1.25 < g < 4 in steps of 0.25. For each {q, L) pair we 
studied all values of k (the number of active colors) in the range 1 < k < [_q\. We thus 
performed 189 runs in total. 

For each triplet (g, fc, L) we performed between 2 x 10^ and 10® CM iterations: more 
precisely, we peformed 10® iterations for 16 < L < 256, 5 x 10'' iterations for L = 512, 
and 2 x 10^ iterations for L = 1024. In terms of the autocorrelation time Tcxp ~ Tint.^-'? 
our total data set at each triplet (g, k, L) ranges in length from 7 x 10''r on the smallest 
lattices at small q (i.e. L = 16, g = 1.25) to ~ 7 x lO^r on the largest lattices at large g 
(i.e. L = 1024, q = 4, k = 1); in 90% of cases (171 out of 189 runs) our run length is at 
least lO^r. These statistics are high enough to permit a high accuracy in our estimates 
of the static (error < 0.1%) and dynamic (error < 1%) quantities, except for the largest 
lattices at large g. 

Our results for the principal static observables, obtained by combining the data for 
all available k values for each pair (g, L), are reported in Tables [3H71 Our results for the 
dynamic quantities Tint,© for the case k = 1 and the most important observables O are 
reported in Tables [HHT^ The complete set of static and dynamic data, for all observables 
and all values of k, is contained as a file cm2data.tar .gz in the preprint version of this 
paper at arXiv.orglll 

The initial configuration of each run was all-bonds-occupied, except for the run at 
{q,k,L) = (3.75,3,1024), for which it was all-bonds- vacant. The number of iterations 
discarded at the beginning of each run in order to allow the system to reach equilibrium 
was j X 10^ for j = 1, 2, 3 or 4, except for eight runs (2 < g < 2.75, L = 256 and A; = 1, 2) 
for which it was 10^ jj Thus, the the discard interval was in all cases at least 100 rint.f and 
in most cases at least 1000 Tint.f- Since our data suggest that r^^p is very close to Tiat,e', 
it follows that our discard interval was in all cases at least 100 Tcxp, which is more than 
sufficient for the systematic error from thermalization to be neglig ibleE 

Our program was written in Fortran 77 and run on a 1266 MHz Pentium III Tualatin 
processor using the g77 Fortran compiler. Our program requires approximately (12 -|- 
26d)L'^ bytes memory for a lattice of linear size L in d dimensions. The CPU time 
required by our program was approximately 0.35 /is/iteration/lattice site; this value is 
weakly dependent on the lattice size and on k and g. In particular, the CPU time per site 
rises sharply on very small lattices due to the "fixed costs" of the algorithm (i.e. those that 

^ The tar file unpacks to make a directory cm2data with subdirectories static and dyncunic. Each 
of these subdirectories is in turn subdivided according to the value of q and, in the dynamic case, the 
value of k. Individual files have names like C1_Q1. 25.txt or tau_Cl_Q4_Kl.txt and have three fields on 
each line: L, value and error bar. 

^ Please don't ask why the discard interval was chosen in this erratic way — we don't remember 
anymore! 

^ Unless there exists a vastly slower mode of which we are unaware, so that in fact Toxp ^ Tint.^'- 
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do not scale with the volume); it rises gradually on lattices L > 128 due to "cache misses" 
coming from the nonlocal nature of the Chayes-Machta algorithm, when the lattice no 
longer fits in the computer's 512 KB cache. The complete set of runs reported in this 
paper used approximately 14 yr CPU time. 

One delicate issue concerns the choice of the pseudo-random- number generator. We 
learned by bitter experience |12] that linear congruential pseudo-random-number gener- 
ators can cause systematic errors in Monte Carlo simulations using the Swendsen-Wang 
(or Chayes-Machta) algorithm, if the lattice size is a multiple of a very large power of 2 
and one random number is used per bond in a periodic manner. These systematic errors 
arise from correlations within a single bond-update half-sweep: see [32] for details. These 
errors can be eliminated (or at least radically reduced) by updating the bonds in a random 
order or in an aperiodic manner, and by using a generator of large modulus (e.g. 60 or 
more bits). In the present project, therefore, we used a linear congruential generator 

Xn+i = aXn + c (mod m) (5.1) 

with modulus m = 2^^, increment c = 1, and multiplier a = 3202034522624059733. This 
multiplier gives good results on the spectral test in low dimensions |13]. Moreover — 
and perhaps even more importantly — we used an "aperiodic" updating, in which the 
random-number subroutine is called only if the two spins are equal and of an "active 
color", i.e. for an edge e = (xy) 

= cTy < k then 
if ran( ) < p then 

else 

He ^ 

endif 
endif 

To the best of our knowledge [12] these choices suffice to make the systematic errors 
negligible. 

6 Data analysis: Static quantities 

Our general methodology for analyzing the Monte Carlo data (both static and dy- 
namic) is as follows: For each quantity Q of interest we impose an ansatz of the form 

Q{L) = AiLP' + A^V' + ... + AkW^ (6.1) 

(or a logarithmic variant thereof), where some of the parameters may be fixed and oth- 
ers free. The precise ansatz will be motivated in each case by finite-size-scaling theory, 
possibly together with some exactly known exponents. Typically our ansatze will have 
2-4 free parameters. We then perform the nonlinear fits corresponding to the chosen 
ansatz by using Mathematica's function NonlinearModelFit. As a precaution against 
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correction-to-scaling terms that we have failed to include in our chosen ansatz, we im- 
pose a lower cutoff L > Lmm on the data points admitted in the fit, and we systematically 
study the the effect on the value of increasing L^i^. In general, our preferred fit for any 
given ansatz corresponds to the smallest L^^^ for which the goodness of fit is reasonable 
and for which subsequent increases in Lmin do not cause the value to drop by vastly 
more than one unit per degree of freedom. In practice, by "reasonable" we mean that the 
confidence level is We do not allow fits with zero degrees of freedom, since 

there would then be no way of testing the goodness of fit. As a last step, we consider 
the effect of including different terms in the ansatz, by comparing the fits obtained with 
different ansatze. 

The exact values of the static critical exponents are of course already known [cf. (11. 5p - 
(ll.lOp ]. Our analysis of the static Monte Carlo data therefore has three major goals: 

1. To test how accurately we would be able to estimate the static critical exponents 
if we did not know them exactly. Here our knowledge of the exact values allows 
us to evaluate systematic errors (including those due to unknown causes) as well as 
statistical ones. This study of the static critical exponents serves as a comparison 
case for our subsequent study of the dynamic critical exponents (Section [7] below). 
In particular, the behavior of the autocorrelation times is closely related to that of 
the specific heat, as is evident from the Li-Sokal bound (11.111) and the empirical 
fact of its near-sharpness. 

2. To estimate universal amplitude ratios such as x* = lim ^{L)/L. 

L—^oo 

3. To study the structure of the corrections to scaling. In particular, are we able to 
see the leading non-analytic correction-to-scaling term L~^^ with exponent given by 
(11.61) ? And what additional corrections are present? Depending on the observable, 
one might expect to see analytic corrections such as L^^ or regular background 
contributions. 

We begin by summarizing the finite-size-scaling (FSS) ansatze that motivate our fits 
(Section 16. ip . We then study successively the correlation length ^ (Section 16. 2p . the 
susceptibility x (Section l6.3p . the largest cluster Ci (Section l6.4p . and the specific heat 
Ch (Section l6.5p . 

6.1 Corrections to finite-size-scaling 

Let us work exactly at the critical point, and consider a quantity Q{L) whose asymp- 
totic behavior as L — )■ oo is 

Q{L) = AL^ (1 + corrections) , (6.2) 



^'■^ "Confidence level" is the probability that would exceed the the observed value, assuming that 
the underlying statistical model is correct. An unusually low confidence level (e.g. less than 5%) thus 
suggests that the underlying statistical model is incorrect] in our context it suggests that the terms we 
have chosen to retain in the ansatz (j6.ip are insufficient to explain our data, i.e. there are unincluded 
corrections to scaling whose contribution is comparable to or larger than the statistical errors in our data. 
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where > is the leading critical exponent and A is the leading critical amplitude 1"1 
Then finite-size-scaling (FSS) theory predicts the following types of contributions to the 
"corrections" in fl6.2p j^1 

1. "Non-analytic" corrections to scaling L^* coming from irrelevant scaling fields, i.e. 
those with i/i < 0. The leading such contribution is given by i/i = yT2 = — •^i 
[cf. fll.6p ]. There are also higher-order "non-analytic" terms of the general form 
]^y,+yj+-^ but these are so far down as to be undetectable in our data. 

2. "Analytic" corrections to scaling coming from nonlinear mixing between the thermal 
and magnetic scaling fields. The leading such contribution behaves as L^^^'^^"^ = 

3. A regular background term [i.e. a contribution of order L° to Q{L)]. 

4. Corrections with negative integer exponents, i.e. L~^, L~^, etc. Some theoretical 
frameworks 01] suggest that these contributions are absent @ 

Looking at the known exact values (ll.5p -( [TTT0|) of these exponents (see Table [2]), we 
can predict the pattern of corrections to scaling for different observables: 

• For strongly divergent quantities such as the susceptibility {ip = •y/u ^ 1.75) and 
the mean cluster size {ip = d — fi/v ~ 1.88), the strongest correction to scaling 
at large L is expected to come from the leading irrelevant term L~^^ whenever 
g > (3 -|- \/5)/2 2.618; indeed, this behavior will hold for all q > 1.294 (resp. 
q > 1.136) if the correction is absent. However, for q < 1.5 the irrelevant, 
regular-background and contributions are all of roughly the same order and 
will be extremely difficult to disentangle. 

• The same pattern is expected to hold for the correlation length, since it is built out 
of the more fundamental quantities x ^"^^ F [cf- fl3.19p ] and therefore inherits their 
corrections to scaling (we expect F to behave like x this context). 

• For weakly divergent quantities such as the specific heat for q > 2, the strongest 
correction to scaling at large L is expected to come from the regular background 
term when 2 < g < 2 + \/2 ^ 3.414 and from the leading irrelevant term when 
g > 2 + ^2. 



For the specific heat we have ip — ajv < when q < 2. The same formulae will apply in this case, 
but the meanings of "dominant" and "subdominant" will be interchanged. 

See e.g. [ini Section 3] for a review and citations to the original literature. We are also indebted to 
Youjin Deng for extremely helpful discussions of these points. 

See [16, Section 3.2] for further discussion. Interestingly, a recent study [45] of percolation on two- 
dimensional lattices, using both transfer- matrix and Monte Carlo methods, found no evidence of any 
correction but found indirect evidence of an correction (namely, an apparent term log L that 
could result from mixing between the and L^^^ = terms when q = 1). 
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• For non- divergent quantities such as the specific heat for g < 2, the dominant 
contribution is the regular background, and the leading singular term LF^I^ becomes 
the principal subdominant contribution. 

One very interesting question concerns the presence or absence of the correction-to- 
scaling term L~^^ = L~^l'^ in the case of the Ising model (g = 2). Plausible physical 
arguments suggest that such a correction should be present in at least some models 
in the Ising universality class, even though special symmetries might cause it to be absent 
in simple exactly-soluble models such as the nearest-neighbor square-lattice model. A 
transfer-matrix study of the square-lattice random-cluster model for q ^ 2 suggested [16] 
that such a correction is indeed present for all g 7^ 2 but has an amplitude that vanishes 
linearly with q — 2 when g — )■ 2. If this is the case, then this correction will be observable 
in derivatives d/dq of standard Potts-model observables evaluated at g = 2. Note also 
that such derivatives can be expressed as bond observables in the random-cluster model; 
therefore, in this scenario, at least some bond observables in the g = 2 random-cluster 
model will exhibit the correction-to-scaling term L~'^^^. 

6.2 The correlation length ^ 

Finite-size-scaling (FSS) theory predicts that C{L)/L for the random-cluster model 
on a torus tends to a universal value x*{q) as L — )■ cxdo More precisely, we expect the 
behavior 

^/L = + BL-P + ... (6.3) 

where p > is one of the correction exponents discussed in the preceding subsection. For 
g = 2, the exact value of x* is known from conformal-invariance theory (17| H8] together 
with numerical integration [39] to be 

x*(g = 2) = 0.9050488292 ± 0.0000000004 . (6.4) 

For g 7^ 2 the exact value of x* is unknown. In this section we will use our Monte Carlo 
data to estimate the universal values x*{q) and to test the agreement of the corrections 
to scaling with the FSS ansatz (16. 3p . 

In Table [15] we report fits of ^/L to the ansatze x* and x* + BL^^ . For g = 1.25 and 
1.50, the fits to both ansatze are excellent, and the estimated values of p agree closely 
with the known exact values of Ai; the correction-to-scaling amplitude B is positive and 
more than three standard deviations away from zero. In particular, there is no evidence 
of any correction. For g = 1.75, the behavior is similar but the quality of the fit is 
poorer (possibly due to a statistical fluctuation that made the L = 1024 point 2a too 
high). For 2 < g < 3.25, the fits to a constant are excellent provided that L^i^ is chosen 
large enough, but the fits to the ansatz x* + BL~^ are mostly badly behaved: either they 
fail completely to converge (as for g = 3.25), or they exhibit estimates for p that seem 
much too low (g = 2) or much too high (g = 2.5,3) or somewhat too high (g = 2.25); 



In general this value depends on the aspect ratio of the torus. Here we are considering only the case 
of aspect ratio 1. 
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only q = 2.75 behaves well in these respects. Also, the estimated correction amplitudes 
B are mostly within two standard deviations of zero. Examination of plots of these fits 
(see Figure d]) suggests that the problem arises because fj L is a very fiat and possibly 
nonomonotonic function of L at large L in this range of This suggests that the ansatz 
(16. 3p might be too simple, and that we are seeing the combined effect of two (or more) 
correction-to-scaling terms of opposite sign. However, our statistical errors are too large 
for us to resolve this question clearly. For 3.5 < < 4, by contrast, the fits to both ansatze 
are again excellent, and the estimated values of p agree well with the known exact values of 
Ai (except at g = 4, where the corrections to scaling log log L/ logL and 1/ logL [12] are 
mimicked by a small inverse power of L). Moreover, the correction-to-scaling amplitude 
B is now clearly negative. This suggests that B passed through zero somewhere in the 
range 2 < g < 3.4 and was small throughout this range, thereby explaining why our fits 
to X* + BL~P were so unstable. 

It is curious that our estimate for x*{q = 2) agrees slightly less well in absolute terms 
with the exact value (16. 4p than the estimate a decade ago by Salas and Sokal [19], despite 
the fact that our raw-data error bars are a factor of 3-7 smaller than theirs. Most likely 
we were the unlucky victim of a statistical fluctuation on our L = 512 and L = 1024 data 
points that placed them roughly 2-3o" and la too low, respectively, thereby causing us to 
choose Lmin = 512 rather than 256 in the fit to a constant and hence to privilege these 
discrepant values [see Figure [Tj^b)]. 

6.3 The susceptibility x 

In Table [TBI we report fits of the susceptibility x to the pure-power-law ansatz AL^^'^, 
just to see how accurately we would be able to estimate the exponent 7/z/ if we did not 
know its exact value. The last two columns show the deviation of the estimated 7/z/ from 
the known exact value (11.101) . in absolute terms and in units of its standard deviation. 

For small q we are able to estimate j/u with extraordinary accuracy: both statistical 
and systematic errors are of order 0.0001 or less. However, as q grows, both the systematic 
and the statistical errors grow: the statistical errors grow because the critical slowing- 
down is becoming gradually worse (see Section [7] below); and the systematic errors grow 
(even faster than the statistical errors) because the corrections to scaling are becoming 
much stronger as the correction-to-scaling exponent Ai decreases (see Table [2]). As a 
consequence, at g ~ 3 the statistical and systematic errors are of order 0.0003 and 0.001, 
respectively; at q = 3.75 the errors are roughly twice this; and at g = 4 the statistical 
and systematic errors are of order 0.0015 and 0.01, respectively. (The large systematic 
errors at g = 4 are of course not surprising in view of the multiplicative logarithmic 
correction (logL)"^/^ and the additive logarithmic corrections loglogL/logL, 1/logL, 
etc. [in]. Indeed, what is somewhat surprising is that the systematic errors are not larger! 
See Section 1^31 below for the case of the specific heat.) 

We next imposed the known exact value (11.101) of •y/u and attempted to extract the 
corrections to scaling. Table [T71 shows fits of x/L^^" to the ansatze A and A + BL~'p. 



Or, perhaps, the statistical errors are just too large, compared to the small "signal", for us to 
estimate the corrections to scaling accurately. 
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Plots for some selected values of q are shown in Figure [2l For q = 1.25 and 1.5, the 
correction to scaling at small L is clearly positive in sign, though very small in absolute 
magnitude (of order 0.001 at L = 16): see Figure [2](a). This suggests that the correction- 
to-scaling exponent is fairly large, in agreement with the theoretical prediction. However, 
the fits produce unusually small estimates for the correction exponent p; this is a possible 
behavior of the "effective exponent" when there are two correction-to-scaling terms of 
opposite sign (e.g. and L~^^). For q = 1.75 and 2, the corrections to scaling are 
small and erratic and it was not possible to fit them to a single inverse power; indeed, 
for q = 2 the behavior may be nonmonotonic in L [see Figure [2|^b)], suggesting again the 
presence of two correction-to-scaling terms of opposite sign. For q > 2.25 the correction 
to scaling is clearly negative in sign, and it gets larger in absolute magnitude as q grows 
(evaluated at L = 16, the magnitude of the correction grows from ^ 0.003 at g = 2.25 to 
~ 0.018 at g = 3.75); furthermore, the fits produce estimates for the correction exponent 
p that are in decent agreement with the exact value of Ai. At g = 4 the correction to 
scaling is extremely strong (as expected) and it is not possible to fit it to an inverse power. 

In summary, for q > 2.25 we obtain modest evidence that a correction L~^^ is present, 
with negative amplitude. For g < 2 we are unable to say much about the corrections to 
scaling except that the amplitude is probably positive for g < 1.5 (we say "probably" 
because it is not clear that the positive contribution at small L comes from the same 
correction-to-scaling term that is dominant at large L). But our data are at least compat- 
ible with the scenario [46J that the correction-to-scaling amplitude passes through zero at 
g = 2. 

6.4 The largest cluster Ci 

In Table [18] we report fits of the mean size of the largest cluster, Ci = (Ci), to the 
pure-power- law ansatz AL'^~^^^, just to see how accurately we would be able to estimate 
the exponent P/u if we did not know its exact value. The last two columns show the 
deviation of the estimated P/u from the known exact value (11.91) . in absolute terms and 
in units of its standard deviation. 

The results are extremely similar to those obtained for the susceptibility. For small 
g, both the statistical and the systematic errors are of order 0.0001 or less. As g grows, 
the statistical errors grow and the systematic errors grow even more. But the results are 
remarkably good except at g = 4. 

6.5 The specific heat Ch 

In Table [19] we report fits of the specific heat Ch to the ansatze AL"^'^ and AL°'^'^ + B. 
The exact value (11.81) of a/u is shown for comparison in the next-to-last column. The 
last two columns show the deviation of the estimated a/u (taken from the fit AL°'^'^ + B) 
from the known exact value, in absolute terms and in units of its standard deviation. 

For g < 2.75 the fits to a pure power law are horrible. This is not surprising: after all, 
for g < 2 the specific heat increases to a finite value as L — )■ oo, so the fit will indicate a 
positive value of a/u (but with a ridiculously poor goodness of fit) while the true value of 
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ajv is actually negative. And for 2 < g < 2.75 the additive constant is very important, 
so a pure power law has poor goodness of fit. 

Good fits to AL'^I'' + B can be obtained for all values of q other than g = 2 by 
appropriate choice of Lmin- For g = 2 the fit fails to converge: initial guesses with A > 0, 
i? < and ajv slightly positive always get driven to a/z^ J, with A ^ —B f +oo. The 
fit is therefore suggesting the correct behavior Cu ~ A\o<gL + B. Indeed, a fit to this 
ansatz is good already for L^^^ = 64 and yields 

A = 0.6371(9), B = 0.1732(42) (x^ = 1.49,3 DF,CL = 68%) , (6.5) 

in good agreement with the known exact values A = 2/7r ^ 0.6366 and B ^ 0.1778 
[SDIEI] If we fit instead to A log^ L + 5 log L + C, we get a good fit for = 32, 

A = -0.0012(6), B = 0.6493(55), C = 0.1415(125) (x^ = 2.36,3 DF, CL = 50%) , 

(6.6) 

in which the extremely small value of A correctly suggests that the log^ L term is absent. 

Although the fits to AL°/^ + B are good for all g 7^ 2, some of the estimated values 
of a/u deviate significantly from the known exact value (see the last two columns of 
Table Hn]). For g < 2 the estimates of a/u are surprisingly good, given that we are 
estimating a subleading singular contribution underneath a nonsingular background. For 
2.25 < g < 2.75 the estimates of a/u are also excellent. However, for g > 3 the estimates 
of a/v deviate by more than three standard deviations from the exact value. This may be 
due to the effect of corrections to scaling: in particular, those governed by the exponent 
Ai, whose known exact value is given by Indeed, for g > 2 + ^2 ^ 3.414 {g > 7/2) 

we have Ai < a/z/, so that the correction-to-scaling contribution is larger than 

the nonsingular background. For 3 < g < 3.75 we therefore tried fits to the ansatze 
AL""/" + SL^/'^-^i and AL""/" + BL"/"-^^ + C, in which Ai is fixed at its known exact 
value and a/u is free: 

(i) For g = 3 it is silly not to include the constant background, since a/u — Ai = 
-0.4 < 0; nevertheless, we obtain a good fit to AL"/" + BL"'/"''^^ when L^in = 64: 

a/u = 0.4152(15), A = 1.397(14), B = -2.346(153) (x^ = 0.75,2 DF, CL = 69%) . 

(6.7) 

However, the estimate for a/z/ is 10 standard deviations away from the correct value! If 
we include the constant background we get a good fit already when L^in = 16, 

a/u = 0.4036(34) A = 1.555(44), B = -0.787(188), C = -0.795(152) 

(X^ = 1.07, 3 DF, CL = 78%) , (6.8) 

and the estimate for a/z/ is in excellent agreement with the correct answer. 

These values of A and B can be extracted from jSU] or [21] after translating their conventions to 

(2) ______ 

ours. We use the definition C\j of specific heat [cf. p.l4p ]. Salas [51 uses the definition C);^' [cf. (15TB ]. 
which is 4 times ours when q = 2. Ferdinand and Fisher use cj^'' multiplied by K^, where K is the 
nearest-neighbor coupling in the Ising normalization; the critical point is given by Kc = 5 log(l + V^). 
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(ii) For q = 3.25 the corrcction-to-scaling contribution is only slightly smaller than the 
nonsingular background, and it will probably be difficult for the fit to separate the two. 
The fit to AL"/^ + BL'^I'"^^ is good already when L^i^ = 16, 

a/v = 0.4972(7) A = 1.201(5) B = -0.945(19) (x^ = 2.58,4 DF,CL = 63%) , 

(6.9) 

but the estimate for a/v is almost six standard deviations away from the correct value 
~ 0.5013. If we include the constant background we also get a good fit when Lmin = 16, 

a/v = 0.4988(44), A = 1.186(42), B = -1.089(396), C = 0.133(365) 

(X^ = 2.45, 3 DF, CL = 49%) , (6.10) 

but the estimate for the constant C is consistent with zero. The estimate for a/u is 
now less than one standard deviation away from the correct value, but this is principally 
because the standard deviation has become much larger, not because the estimated value 
has actually moved much closer to the true value! 

(iii) For q = 3.5 the correction-to-scaling contribution is now slightly larger than the 
nonsingular background, and it will again likely be difficult for the fit to separate the two. 
The fit to AL""/" + BL"^'"-^^ is good when L^^^ = 64, 

a/u = 0.5945(39), A = 0.926(25), B = -0.017(91) (x^ = 0.91,2 DF,CL = 64%) , 

(6.11) 

but the estimate for the correction-to-scaling amplitude B is consistent with zero. Indeed, 

the result from this fit is virtually identical to what was obtained from the ansatz AL"^'^ + 

B, and the estimate for a/u is again about four standard deviations away from the correct 
value ^ 0.6101. By contrast, the fit to AL'^/" + 5L"/^-^i + C is good already when 

Lmin 16, 

a/u = 0.6095(67), A = 0.805(43), B = 2.961(368), C = -3.509(390) 

(X^ = 2.32, 3 DF, CL = 51%) , (6.12) 

and the estimate for a/u is now in excellent agreement with the correct answer. Inter- 
estingly, the estimates for the amplitudes B and C are not consistent with zero; rather, 
they are strongly nonzero but of opposite signs. Clearly, what happened is that when we 
performed a fit with a single correction term (whether B or BL"^'^~'^^ , which are anyway 
nearly the same) this pair of correction terms combined to make an "effective" correction 
term (in the given range of L) with a nearly zero amplitude; but this gave a biased esti- 
mate of the leading exponent a/z/. What is slightly surprising is that our fit was able to 
separate the correction-to-scaling contribution iJ/,~0 090 f^Qj^i the nonsingular background 

C. Perhaps the stunning agreement of the estimated value oi a/p with the exact answer 
is a fluke and ought not be taken too seriously. 

(iv) For q — 3.75 the correction-to-scaling contribution is significantly larger than the 
nonsingular background. The fit to AL°'/'^ + BL°'/^~^'^ is good when Lmin = 64, 

a/u = 0.7181(90), A = 0.568(42), B = 0.413(76) (x^ = 1.26,2 DF,CL = 53%) , 

(6.13) 
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and the estimate for a/z/ is only about two standard deviations away from the true value 
^ 0.7376. By contrast, the fit to AL""/" + 31°"/"-^^ + C is good already when L^-,^ = 16: 

a/u = 0.7856(237), A = 0.288(68), B = 1.085(64), C = -1.046(37) 

(X^ = 0.30, 3 DF, CL = 96%) . (6.14) 

Here the estimated amplitudes B and C have opposite signs and are apparently nonzero; 
but the estimate for a/u has now far overshot the correct value (it is again two standard 
deviations away, but with a much larger standard deviation). The poor performance of 
this two-correction-term fit — in a case where the two terms BL~^'^^^ and C should have 
been much easier to separate than they were for q = 3.5 — suggests that the good result 
obtained for q = 3.5 was indeed a fluke and that the mediocre result obtained for q = 3.75 
is what should ordinarily be expected. 

Finally, for g = 4 the true leading behavior is known [T9] to be L(logL)~^/^, but with 
corrections to scaling down by log log L/ log L, 1/logL, etc. It is clearly hopeless to try 
to fit to such an ansatz unless one has data for colossally large values of L. It is of course 
not surprising that fits to AL°'^'^ or AL°'^'^ + B gave estimates of a/u near 0.8, far off from 
the correct value 1; the factor (logL)'^/^ is imitating a power in our range of L. 

7 Data analysis: Dynamic quantities 

In this section we analyze the dynamic data by the same general methods as were 
used in the preceding section to analyze the static data. Our main goal is to estimate the 
dynamic critical exponents Zmt.o associated to the integrated autocorrelation times Ti^t,o 
for various observables O. 

We proceed as follows: First we discuss the dependence of the autocorrelation times 
on the number k of active colors (Section l7.ip . and we give an overview of the qualitative 
behavior of the autocorrelation times for different observables (Section 17.21) . Then we 
present a detailed analysis of the dynamic critical exponent Zij^t,e' (Section [7.31) : in partic- 
ular we discuss the sharpness of the Li-Sokal bound (Section I7.4p and the correctness of 
the Ossola-Sokal conjecture (Section 17.51) . Finally, we analyze briefly the dynamic critical 
exponent -Zint,© for other observables O (Section l7.6p . 

7.1 Dependence on k 

We began by analyzing the dependence of Ti^t,o on the number k of active colors. Of 
course we expect that all values of k lie in the same dynamic universality class: that is, we 
expect that the ratios of Tint,© for different k tend to nonzero finite constants as L — )■ cxd. 
Moreover, it is intuitively reasonable to think that an update with k active colors does 
roughly "fc times as much work" as an update with one active color; therefore, we expect 
that Tint,o should be roughly proportional to 1/k. 

We tested these expectations by analyzing the ratios Tint,c'(l)/Tint,c'(fc) as a function 
of L for each (g, O) and each allowable k. In all cases the ratios are fairly close to k, as 
expected; but in general they are not exactly equal to from k. The ratios also show some 
dependence on L, but tend as L — )■ oo to a limiting value, again as expected. Roughly 
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speaking, for the smaller values of q the L-dependence is fairly strong, and the ratios are 
comparatively far from k (which is perhaps not surprising because the values of Tint,c' cire 
themselves quite small); the limiting values also appear to be different from k, though this 
conclusion is only tentative because of the strong corrections to scaling. For the larger 
values of q the L-dependence is weaker, and the ratios are closer to fc; in particular, the 
limiting values are compatible with k within our statistical errors. In Table [20] we show 
typical examples of these two behaviors, namely O = £' with q = 2 and q = 3.25. In 
Table 121] we show our best estimates for the limiting ratios T"int,£-'(1)/Tint,f'(fc), obtained 
by fitting the ratio to a constant and increasing Lmin until a decent fit is obtained. The 
behavior for the other observables is qualitatively similar. 

Having confirmed that all values of k lie in the same dynamic universality class, we 
henceforth analyze the data for each value of k separately and then compute a weighted 
average of the resulting exponent estimates. 

7.2 Summary of qualitative behavior 

Let us begin by summarizing the qualitative behavior of Tint,© for different observables 

O: 

1) For nearly every triplet (g, k, L), we find that S' is the observable (of those we have 
measured) that has the largest Tint- The only exceptions are g = 4, L = 16, 1 < A; < 4 
(for which Xint.Sa is slightly larger than Ti^t,e') and q = 4, L = 1024, k = 1 (for which 
Tint.Ci is slightly larger than rint,£-')- But these differences are extremely small and may 
well represent statistical fluctuations. 

2) For every triplet {q,k,L), we find that C2 is the observable (of those we have 
measured) that has the smallest Tint- 

3) For every triplet (g, k, L) we find that 

Tint.f > Tint,Ar • (7.1) 

When q is an integer and k = q, this inequality is easily proved rigorously for the SW 
algorithm |T7], so it is not surprising that it holds here for CM. However, we do not yet 
have a rigorous proof (not even in the case when q is an integer and k < q). 

4) In general, the observables O that we have measured fall into four groups according 
to their integrated autocorrelation times Ti^t,o- 

1. S', A/", ^2 and Ci have the largest values of Tint, and they are all fairly close to each 
other (all are at least ^ 0.75 times that of S', and usually much closer); 

2. J-"' has an intermediate value of Tint, of order 0.5-0.9 times that of S'; 

3. C3 has a slightly lower value of Tint, of order 0.4-0.8 times that of S'; 

4. C2 has the smallest Tint, of order 0.35-0.8 times that of S'. 

Indeed, for all triplets (g, k, L) we have Ti^t,Oi > Tint,c»2 > ''"mt.Os > Tint,04 whenever 
Ci, (92,(^3,(94 belong to groups 1,2,3,4, with the exception that for a few triplets at 
g = 4, L > 256 we have Ti^t,T' < ^int.Cs- 
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These behaviors can be better understood by looking at the normahzed autocorrelation 
functions pooif). A typical example is shown in Figure |3l We see that pe'e'it) is nearly 
a pure exponential, so that Tij^t,e' ~ ^exp- By contrast, the autocorrelation functions for 
the observables O in groups 2, 3 and 4 exhibit an initial fast decay, followed by a decay 
at the same exponential rate Texp but with an amplitude Ao that is significantly less 
than 1. What we do not know is whether Ap tends to a nonzero value as L — )■ oo (in 
which case we will have Zi^t^o = ^exp) or tends to zero as an inverse power of L (in which 
case we will have 2int,o < -^exp)- See Section for further analysis of this question; 
and see [201 Section 5.2] for a more detailed analysis in the case of the Swendsen-Wang 
dynamics for the three-dimensional Ising model. 

In the following subsections we shall fit Tint,©, for each observable O, to a variety of 
ansatze, notably: 

• Fits for z = 0: nr,t,o = A or Ti^t,o = A + BL~p. 

• Fits for z = with a multiplicative logarithm: Tiat,o = ^logL + B or rint,o = 
A\og^ L + BlogL + C. 

• Fits for z > 0: Tint,o = oi Tint.o = + B. 

Note that the fits to A + BL~^ and AL^ + B are in fact the same fit in different notation! 

We shall begin (Sections 17.3117.51) by focusing on the observable £', which has the 
largest autocorrelation time of all the observables we measured. Then (Section I7.6p we 
shall discuss, more briefly, the other observables, with emphasis on C2. 

7.3 Dynamic critical exponent Zint,£' 

In this subsection we flt the integrated autocorrelation time T[^t,e' to a variety of ansatze 
in an effort to estimate the dynamic critical exponent -Zmt.f- We begin by presenting our 
flts, in order of increasing q, without much comment. Then we go back and try to 
interpret what these flts might be telling us about the dynamic critical behavior of the 
Chayes-Machta algorithm as a function of q. 

For q = 1.25 the behavior is fairly clear: Tij^t,s' tends to a flnite constant as L — )■ 00. 
The flts to the ansatz Tint,£' = A are horrible (except of course the flt Lmin = 1024 that 
has zero degrees of freedom); but if we flt to Ti^rt,£' = A + BL~p we get a decent flt for 
Lmin = 128: 

A = 2.19(8), B = -1.61(5), p = 0.213(35) (x^ = 0.49, 1 DF, CL = 48%) . (7.2) 

Finally, the flt to AL^ yields an estimated exponent z ~ 0.05 but has poor goodness of 
flt even when L^in = 256: z = 0.053(1), A = 1.27(1) with = 3.76, 1 DF, CL = 5%. 

For q = 1.5 the flts to a constant A are again horrible; but if we flt to A + BL^^ we 
get a good flt already for Lmin = 32: 

A = 17.09(2.96), B = -15.86(2.92), p = 0.034(7) (x^ = 0.82,3 DF,CL = 85%) . 

(7.3) 
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We also tried a fit to AlogL + B: a. decent fit is obtained for Lmin = 128, namely 

A = 0.440(4), 5 = 1.486(21) (x^ = 0.83, 2 DF, CL = 66%) . (7.4) 

A good fit to A\og^ L + BlogL + C is obtained already for L^in = 32, but with a 
slightly negative value of A (which is of course impossible for the actual asymptotics) : 
A = -0.008(2), B = 0.53(2), C = 1.24(4) with = 0.80, 3 DF, CL = 85%. Finally, the 
fit to AL^ yields an estimated exponent z ~ 0.10 but again has poor goodness of fit even 
when Lmin = 256: z = 0.104(2), A = 2.20(2) with = 3.84, 1 DF, CL = 5%. 

For q = 1.75 (and all larger q) the fit to a constant A is again horrible. The fit to 
AhgL + B is mediocre even when L^i^ = 256: A = 1.36(2), B = 0.32(12) with = 2.73, 
1 DF, CL = 10%. A decent fit to A log^ L + B\ogL + Cis obtained already for L^in = 16: 
A = 0.037(3), B = 0.915(23), C = 1.689(45) with = 4.14, 4 DF, CL = 39%. The fit 
to AL' is good when L^in = 256: z = 0.158(2), A = 3.28(5) with = 0.21, 1 DF, CL = 
65%. The fit to AL' + B is decent already for L^m = 16 (x^ = 4.63, 4 DF, CL = 33%), 
but the x^ drops notably from Lmin = 64 (x^ = 2.86) to Lmin = 128 (x^ = 0.33), so that 
our preferred fit is = 128: 

z = 0.085(26), A = 9.66(4.42), B = -7.56(4.82) (x^ = 0.33,1 DF, CL = 57%) . 

(7.5) 

For g > 2 we have data from more than one value of k. For simplicity we discuss here 
in words only the case k = 1; the other cases are qualitatively similar and are reported in 
Tables [22] and El 

For g = 2 the fit to AlogL + 5 is poor (CL < 2%) even when L^in = 256. The fit 
to Ahg^ L + BhgL + C is good already for L^.^ = IQ: A = 0.206(6), B = 1.046(53), 
C = 2.572(103) with x^ = 1-40, 4 DF, CL = 84%. The fit to AL' is good for L^in = 256: 
z = 0.215(3), A = 4.47(9) with x^ = 0.06, 1 DF, CL = 81%. Finally, the fit to AL' + B 
is good already for Lmin = 16: 

z = 0.145(4), A = 10.4(5), B = -8.4(6) (x^ = 2.32,4 DF,CL = 68%) . (7.6) 

For q = 2.25 the fit to A log L + B is again poor (CL < 5%) even when Lmin = 256. The 
fit to Alog^L + SlogL + C is decent already for L^^in = 32: A = 0.77(3), B = -0.84(24), 
C = 7.07(53) with x^ = 3.15, 3 DF, CL = 37%. The fit to AL' yields an estimated 
exponent z ~ 0.29 but has mediocre goodness of fit even when L^i^ = 256: z = 0.286(4), 
A = 5.35(14) with x^ = 2.16, 1 DF, CL = 14%. Finally, the fit to AL^ + B is good already 

for -^rnin ~ 32 i 

z = 0.235(8), A = 8.91(61), B = -6.68(86) (x^ = 1.62,3 DF,CL = 66%) . 

(7.7) 

For q = 2.5 (and all larger q) the fit to AlogL + i? is poor (CL < 3%) even when 
L^in = 256. The fit to Alog^L + 51ogL + C is decent for L^.^ = 128: A = 2.7(3), 
B = -13.7(3.1), C = 37.4(8.6) with x^ = 0.66, 1 DF, CL = 42%. However, the error 
bars on B and C are very large in absolute magnitude (and in particular large compared 
to the value of A), which renders the fit somewhat dubious. The fit to AL^ is good for 
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Lmin = 256: z = 0.353(6), A = 6.35(23) with = 0.02, 1 DF, CL = 90%. Finally, the fit 
to AL^ + B is excellent already for Lmin = 32: 

z = 0.315(8), A = 9.01(59), B = -6.75(1.00) (x^ = 0.21, 3 DF, CL = 98%) . 

(7.8) 

For q = 2.75 the fit to Alog^ L + BlogL + C is decent for L„un = 64: A = 5.9(3), 
B = -31.9(2.7), C = 71.0(6.5) with = 2.56, 2 DF, CL = 28%. But the error bars 
on B and C are again quite large. The fit to AL^ is good for L^nm — 256: z — 0.424(8), 
A = 7.19(34) with = 0.07, 1 DF, CL = 80%. Finally, the fit to AL^ + B is excellent 
already for L^in — 16: 

z = 0.411(5), A = 8.19(30), B = -4.76(57) (x^ = 1.06,4 DF,CL = 90%) . 

(7.9) 

For q — 3 the fit to A log^ L + B log L + C is decent for Lj^m = 128, but with huge 
error bars: A = 16.0(1.5), B = -112.6(16.5), C = 257.5(44.5) with x^ = 0.77, 1 DF, CL 
= 38%. (The same behavior persists for all larger values of q. with even larger values of 
the coefficients A,B,C and their error bars; we refrain from reporting the gory results.) 
The fit to AL^ is good for L^nm = 128: z = 0.505(6), A = 7.52(23) with x^ = 1-18, 2 DF, 
CL = 56%. Finally, the fit to AL^ + B is excellent already for Lmm = 32: 

z = 0.481(10), A = 9.08(63), B = -6.79(1.68) (x^ = 0.84,3 DF,CL = 84%) . 

(7.10) 

For q = 3.25 the fit to AL^ is decent for L^in = 128: z = 0.590(7), A = 7.66(30) with 
X^ = 3.14, 2 DF, CL = 21%. Finally, the fit to AL' + B is decent for L^i^ = 64: 

z = 0.558(21), A = 9.8(1.4), B = -11.9(6.0) (x^ = 2.93,2 DF,CL = 23%) . 

(7.11) 

For q = 3.5 the fit to AL' is decent for L^^^ = 128: z = 0.676(10), A = 7.77(39) with 
X^ = 1.07, 2 DF, CL = 58%. Finally, the fit to AL^ + B is decent for L^in = 64: 

z = 0.648(25), A = 9.5(1.6), B = -14.4(8.9) (x^ = 1.15,2 DF,CL = 56%) . 

(7.12) 

For q = 3.75 the fit to AU is good already for L„,i,, = 32: z = 0.779(4), A = 7.14(13) 
with x^ = 0.94, 4 DF, CL = 92%. Finally, the fit to AL' + B is excellent already for 

^min = 16: 

z = 0.790(9), A = 6.7(3), B = 3.0(1.5) (x^ = 0.33,4 DF,CL = 99%) . (7.13) 

For g = 4 the fit to AL^ is good for L^in = 128: z = 0.916(16), A = 5.44(46) with 
= 0.15, 2 DF, CL = 93%. Finally, the fit to AL^ + B is excellent already for L^^in = 32: 

z = 0.935(18), A = 4.8(5), B = 20.5(5.1) (x^ = 0.66,3 DF, CL = 88%) . (7.14) 

Let us now comment on what we think these fits show. 

For q = 1.25 it seems fairly clear that rint.f converges as L — )■ oo to a finite value, i.e. 
Zmt,e' = 0. For q = 1.5 the behavior is unclear: perhaps Tint.f converges to a finite value, 
but with extremely strong corrections to scaling (e.g. A + BL~p with p > extremely 
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small); or perhaps it grows logarithmically, like logL or log^ L; or perhaps it grows with 
a very small positive power {z < 0.1) together with very strong corrections to scaling. 
The evidence points weakly towards the first scenario, but a divergence like logL or an 
extremely small positive power of L is also a possibility. For 1.75 < q ^ 2.25 a reasonable 
fit is obtained with the ansatz A\og^ L + B logL + C; but it seems to us implausible on 
theoretical grounds that we would have such a logarithmic growth for an entire interval 
of q. Much more likely is that there exists one value such that z = for q < q^ and 
2; > for g > g^, in which case Ti-at,e' might grow logarithmically at q = q^, (but only 
there). Our data suggest that q^ lies between 1.25 and 1.75; our best guess would be 
~ 1.6, based on linearly interpolating the exponent estimates produced by our fits for 
q = 1.5 and q = 1.75. 

For all g > 1.75 we are able to obtain decent fits to the ansatze AL^ and AL^ + B 
with an exponent 2; > 0; these fits are reported in Tables [22] and [231 respectively, and 
the results obtained by averaging over k are reported in Table [231 For the smaller values 
of g in this table, the discrepancy between the fits with and without a constant term B 
is fairly large: this is not surprising because z is fairly small and hence the effect of the 
constant term is very strong. As g grows, the discrepancy between the two fits decreases: 
from ^ 0.1 at g = 1.75 and ~ 0.07 at g = 2 to approximately zero at g = 3.5; for g > 3.5 
the discrepancy has the opposite sign but remains small. Correspondingly, the estimated 
value of B appears to go through zero (and change sign) at g ~ 3.5. We have a slight 
preference for the fits to AL^ + B, for the simple reason that such a constant term must 
surely be present, if only because the definition of Tint,© is somewhat arbitrary (should one 
include the contribution | from t = or not?). But our data are insufficient to resolve 
clearly the discrepancy between the two fits. We therefore choose to report our results 
for g > 1.75 in the form 

best estimate ± statistical error ± systematic error , (7-15) 

where "statistical error" denotes the one-standard-deviation error bar from the fit to 
AL^ + B (after averaging over k the fits from the chosen values of Lmin); and "systematic 
error" is a 68% subjective confidence interval defined as the absolute value of the discrep- 
ancy between the fits AL^ and AL^ + B (after averaging over k the fits from the chosen 
values of L^m) plus 0.02. The final results are therefore: 



q = 


1.75 




= 0.085 ±0.026 ±0.094 


(7.16a) 


q = 


2.00 


Zmt,S' 


= 0.143 ±0.003 ±0.092 


(7.16b) 


q = 


2.25 


Zmt,S' 


= 0.231 ±0.008 ±0.071 


(7.16c) 


q = 


2.50 


Zint,e' 


= 0.307 ±0.007 ±0.066 


(7.16d) 


q = 


2.75 


Z\-at.,e' 


= 0.408 ± 0.005 ± 0.036 


(7.16e) 


q = 


3.00 




= 0.497 ±0.003 ±0.033 


(7.16f) 


q = 


3.25 


Zint,£' 


= 0.572 ±0.007 ±0.035 


(7.16g) 


q = 


3.50 


Zint,e' 


= 0.689 ±0.004 ±0.024 


(7.16h) 


q = 


3.75 


Zint,e' 


= 0.796 ±0.004 ±0.032 


(7.16i) 


q = 


4.00 


Zint,S' 


= 0.910 ±0.005 ±0.032 


(r.i6j) 
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In Figure m we plot these estimates versus q, and compare them with the static exponents 
a/u and P/u. The Ossola-Sokal conjecture z> 13/v appears to be violated for 1 < g < 
1.95 (but see Section [731 for an alternative fit that is compatible with the conjecture). 
The Li-Sokal bound z > a/z/ is obeyed for g 7^ 4 and appears to be non-sharp (see 
Section [73] for a more detailed analysis). The apparent violation of the Li-Sokal bound 
at g = 4 is manifestly due to the multiplicative logarithmic corrections: it is known [19] 
that Ch ~ L(logL)-3/2 ^^^^ 

a power law (with or without a constant background) 
yield an effective exponent ~ 0.8 (see the last line of Table [TIJl) : it is therefore not surprising 
that Tint,£' shows a similar behavior. If one looks directly at the ratio Tint,£'/^-ff finds 
that the Li-Sokal bound (which is after all a rigorous theorem!) is obeyed (see Section [73]). 

For q = 1.25 and 1.50 our best estimates suggest that n^t^s' is bounded as L — i- 00, 
i.e. that Zint,£' = 0. Our fits to A + BL~p suggest values for the correction exponent p [cf. 
f l7.2p / ( TTS]) ]. but we do not know how reliable these estimates are. 

Since the autocorrelation function ps'£'{t) is nearly a pure exponential, we expect that 
the dynamic critical exponent Zexp is either exactly equal or almost exactly equal to Zint,£'- 

7.4 Sharpness of Li— Sokal bound 

The estimates of Zint.f summarized in Figure [1] suggest that the Li-Sokal bound 
Zmt,£' > a/^ holds as a strict inequality over the entire range 1 < g < 4, i.e. that it 
is non-sharp by a power. But this conclusion is weakened by the fact that our fits of the 
specific heat Ch give estimates of a/u that deviate significantly from the known exact 
values (Table UM. It is therefore of interest to study directly the ratio Tint,£'/CH, in an 
effort to fit its behavior as L — )■ 00 to one of the following ansatze: 

1) Asymptotically constant with additive corrections to scaling A + BL~^. 

2) A logarithmic growth, either as A\ogL + S or as A\o^ L. 

3) A power-law growth AL^ or AL^ + B. 

Unfortunately our time-series analysis does not produce statistically valid error bars 
for composite static-dynamic quantities such as Ti^t,£'/CH- We therefore conducted the 
analysis in this subsection under the crude assumption that the statistical fluctuations on 
our estimators of rint,£-' and Ch are uncorrelated. In fact it is likely that these fluctuations 
are positively correlated, so that the true error bars on the ratio Ti^^t,£'/CH are smaller 
than we have supposed. If so, this means that the true confidence level of our fits is 
smaller than what we report]^ 

The fits to A + BL~'p with p > were always bad: either they had a horrible confidence 
level, or they converged to a value p < (indicating that the true leading behavior is 
a power- law growth). We therefore focussed on comparing the logarithmic-growth and 
power-law-growth scenarios. We chose AlogL + B and AL'^ as the ansatze in order 
to compare fits with an equal number of free parameters. The results are shown in 

This failure to provide statistically valid error bars is embarrassing. We should have used the 
batch-means method 1 171 Section 4.2] to obtain such error bars. Unfortunately, the raw data from our 
simulations are no longer accessible, so we are unable to conduct such a reanalysis. 
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Table [251 We see that both fits are in general good (though the confidence levels may 
be overestimated as noted above). However, the power-law fits are in general better. In 
particular, there are no cases in which the power-law fit has a confidence level less than 
25%; but there are quite a few cases in which the logarithmic fit exhibits such a low 
confidence level (sometimes much lower). We therefore conclude that the fits to Ti^t,£' /Ch 
also provide weak evidence that the Li-Sokal bound is non-sharp by a power. 

7.5 Test of Ossola— Sokal conjecture 

In the Introduction we argued that the Ossola-Sokal conjecture zqm > /3/z/ is probably 
false, on the grounds that it fails for g = 1 (where zcm = and P/u > and that zcm is 
presumably a continuous function of q. But this latter assumption is far from certain: it 
is possible, for instance, that 2cm = exactly for all q near 1, but with an amplitude 
that vanishes as g J, 1. In this subsection we would like to test this scenario against our 
data for 1.25 < q < 2. 

We first tried fits to rj^t.^' = AL^/^+B. We then tried fits to Tint,£' = AL'^/'^+B+CL-p 
where p > is fixed and A, B, C are free. Let us report our results from these fits in 
decreasing order of q: 

1) For g = 2 (A; = 1) the fit to ALf^/" + B is decent already for L^^^ = 32 (x^ = 3.90, 
4 DF, CL = 42%) and better for L^nm = 64: 

A = 13.27(7), 5 = -11.83(13) (x^ = 1.51, 3 DF, CL = 68%) . (7.17) 

This behavior is not surprising, because our fit to AL^ + B with z free yielded z = 0.145(4) 
with Linin = 16, which is not very far from (S/u = 0.125. 

In Figure [5] we plot the results of our fits to AL^^'^ + B + CL~'p as a function of p, for 
Lmin = 16 and 32: 

• For Ljnin = 16 an excellent confidence level is obtained over the whole range < 
p < 2 (and indeed beyond), with an optimum at p ~ 0.692 (CL = 92%). The 
estimated amplitudes A are positive and very far from zero; at the optimum we 
have A = 13.44(9). As p — )■ oo the amplitude A tends to a value ~ 13, which is 
close to that obtained from the fit to AL^^''' + B. 

• For Lmin = 32 an excellent confidence level is again obtained over the whole range 
< p < 2 (and indeed beyond), with an optimum at p ~ 0.228 (CL = 84%) but 
with a very broad peak. The estimated amplitudes A are again positive and very 
far from zero; at the optimum we have A = 13.71(30). As p — t- oo the amplitude A 
again tends to ~ 13. 

2) For q = 1.75 the fit to AL^/^ + B is decent already for = 128 (x^ = 2.31, 2 DF, 
CL = 32%) and excellent for L^in = 256: 

A = 5.36(8), 5 = -2.62(17) (x^ = 0.002, 1 DF, CL = 97%) . (7.18) 

Once again this behavior is not surprising, because our fit to AL^ + B with z free yielded 
z = 0.085(26) with L^in = 128, which is not very far from (3/u ^ 0.121. 
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In Figure [6] we plot the results of our fits to AL^^'^ + B + CL ^ as a function of p, for 
-^min = 32 and 64: 

• For Lmin = 32 a decent confidence level is obtained for p < 0.6, with an optimum 
at p ~ 0.056 (CL = 41%). The estimated amplitudes A are positive and far from 
zero; at the optimum we have A = 3.51(24). As p — )■ oo the amplitude A tends to 
a value ~ 5.8, which is close to that obtained from the fit to AL^^'^ + B. 

• For Lmin = 64 a decent confidence level is obtained over the entire range < p < 2 
(or even beyond), with an optimum at p ~ 1.321 (CL = 96%). The estimated 
amplitudes A are again positive and far from zero, and in fact quite close to those 
obtained from Lmin = 32. At the optimum, we have A = 5.28(7). 

3) For q = 1.5 the fit to ALl^'^ + B is poor even for L^^^ = 256 (CL = 3%). In FigureE] 
we plot the results of our fits to AL^^'^ + B + CL~^ as a function of p, for Lmin = 32 and 
64: 

• For Ljnin = 32 a decent confidence level is obtained for p < 0.3, but with the 
optimum attained at the ridiculously small value p ~ 0.004 (CL = 86%). Also, 
the estimated amplitude A is negative for p < 0.033, which is obviously impossible. 
However, for 0.08 ^ p ^ 0.3 a decent fit is obtained with an amplitude A that is 
positive and far from zero. 

• For Ljnin = 64 a decent confidence level is obtained for < p < 0.9, with an 
optimum at p ^ 0.228 (CL = 86%). The estimated amplitudes A are again negative 
for p < 0.038 but positive and far from zero for p > 0.1. At the optimum, we have 
A = 1.11(8). 

4) For q = 1.25 the fit to AL^/'' + B is poor even for L^ain = 256 (CL = 1%). The 
fits to AL'^/" + B + CL-P are poor (CL < 6%) for L^,^ = 32 and 64 for all p > 0. For 
Lmin = 128, however, we are able to obtain decent fits over the whole range < p < 2, 
as shown in Figure [HI The optimum lies at p ~ 0.695 (CL = 99.9997%). The estimated 
amplitudes A are negative for p < 0.211 but positive and far from zero for p > 0.4. At 
the optimum, we have A = 0.28(2). 

The foregoing fits show that the Ossola-Sokal conjecture zqm > fS/u is not ruled 
out by our data at 1.25 < q < 2. Indeed, our data are consistent with the possibility 
that zqm = P/^ exactly for 1 < g < 2, but with an amplitude that vanishes as g J, 1, 
perhaps proportional to g — 1 (the fits for q = 1.25 and 1.50 are consistent with this latter 
behavior). 

7.6 Dynamic critical exponents Zint,o for other O 

Let us now look briefly at the dynamic critical exponents -Zint.o for observables O other 
than S'. For O = J\f,S2,Ci, the values of Tint,© are very close to those of nnt^s', so the 
estimates of ^mt,© will be nearly the same; little would be gained by going through these 
fits in detail. Instead, it seems sensible to look at the observable that has the smallest 
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autocorrelation time, namely C2 — the idea being that if any differences in Zmt.o between 
different observables are to be detected, they will be detected here. 

From Figure [3] we see that the autocorrelation function of S' is nearly a pure exponen- 
tial (and this is so for all q, k, L), so that Ti^t,s' ~ ^exp and hence -Zmt.f = -^exp or nearly 
so. The autocorrelation function of C2, by contrast, exhibits an initial fast decay, followed 
by a decay at the same exponential rate Texp as for the other observables but with an 
amplitude Ac^ that is significantly less than 1 (e.g. around 0.4 in the plot shown). The 
key question is whether this amplitude tends to a nonzero value when L — )■ cxd (in which 
case we will have Zint,C2 = Aixp) or tends to zero as an inverse power of L (in which case 



In Tables [26] and [27| we show the fits for Ti^t,C2 ^o the ansatze AL^ and AL^ + B. In 
Table [28] we show the results for Zint,C2 obtained by averaging over k. 

The estimates for Zint,C2 are indeed less than those for Zmt.f, by an amount that is 
fa 0.1 for q ^ 2 and tends to zero as g ^ 4. The question is: Are these differences real, 
or are they artifacts of corrections to scaling at small L? The fact that the differences 
are smaller for g > 3, where the autocorrelation times are larger, suggests that perhaps 
the differences will disappear as L — )■ 00 but that at small q we have to go to larger L to 
see this. But this is far from clear; we will only know the truth by doing simulations at 
significantly larger values of L. 

8 Discussion 

In this paper we have studied the dynamic critical behavior of the Chayes-Machta 
algorithm as a function of q over the whole range 1.25 < q < 4. We have obtained 
estimates of the dynamic critical exponent -Zmt.f as a function of q: see fl7.16p . Table 12^ 
and Figure H] Since the autocorrelation function p£i£i (t) is nearly a pure exponential, 
we also expect that the dynamic critical exponent Zexp is either exactly equal or almost 
exactly equal to Zi^t,e'- 

By simultaneously studying the whole range of values of q, we were able to gain some 
insights that would not have been available had we studied only a single value of q (such 
as the Ising value g = 2) or even all integer values of q as in past studies of the Swendsen- 
Wang algorithm. For instance: 

1) The autocorrelation time Ti^t^e' at g = 2 can be plausibly fit with the ansatz 
A log^ L + B log L + C, suggesting that the Li-Sokal bound might be sharp modulo a 
logarithm; this agrees with the conclusions of the paper [16], where the data were in 
fact found to slightly favor the non-sharp-by-a-logarithm ansatz over the non-sharp-by-a- 
power ansatz. But we now find that the good fit to A log^ L + B log L + C persists over the 
whole range 1.75 < g ^ 2.25. And it seems to us implausible on theoretical grounds that 
we would have such a logarithmic growth for an entire interval of g; rather, we expect 
that there exists one value g^ such that z = for q < q^ and z > for q > q^, with a 
possible (poly) logarithmic growth at q = q^. Our data suggest that g* lies between 1.25 



^® See Section 5.2] for a more detailed analysis of this kind, for the Swendsen-Wang dynamics for 
the three-dimensional Ising model. 
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and 1.75, with our best guess being around 1.6. If this scenario is correct, it follows that 
the Li-Sokal bound is non-sharp by a power for < g < 2. This then suggests (but 
of course does not prove) that it might be non-sharp by a power over the whole interval 

< g < 4 (possibly reverting to non-sharp by a logarithm at q = 4). 

2) By considering the Ossola-Sokal conjecture ^cm > simultaneously for all g, we 
can see immediately that it fails at g = 1 (where zqm = but /3/z/ > 0) and hence fails also 
for g near 1 if the dynamic critical exponent zcm is a continuous function of g. Indeed, 
our pure power-law fits suggest that the conjecture fails for 1 < g < 1.95. (In particular, if 
g^, > 1 as just suggested, then the conjecture fails spectacularly in the interval 1 < g < g^,, 
where zqm = but /3/z/ > 0.) However — and perhaps surprisingly — our data are also 
compatible with an alternative scenario in which zqm = exactly for all g near 1, but 
with an amplitude that vanishes as g J, 1 (Section I7.5p . 

However, the behavior of the Chayes-Machta autocorrelation time for 1 < g < 2 
is still unclear: though our data strongly suggest that rint,£-' diverges as L — )• oo for g 
slightly below 2, and is nondivergent as L — )■ oo for g slightly above 1, we cannot rule out 
the possibility that the data at L < 1024 are misleading and that the true asymptotic 
behavior is different from what we conjecture. Future work at larger values of L would of 
course be desirable. 

Likewise, though our data suggest that the Li-Sokal bound is non-sharp by a power 
for 1.6 < g < 4, it is also true that the exponent estimates have been dropping over time 
as data becomes available at larger values of L and as we try ansatze other than a pure 
power law (see Table |2^ for the effect of the ansatz). Our data suggest most clearly the 
non-sharpness of the Li-Sokal bound (when the non-sharpness is measured in units of 
the standard deviation of our estimate) at g = 2.75 and g = 3. At these values of g, 
the critical slowing-down is strong enough that interference from the regular background 
term is less important than it is at smaller g (i.e. the AL^ and AL^ + B fits show less 
discrepancy) but modest enough that we can have reasonably good data on large lattices 
(contrary to the situation at larger g); furthermore, the correction-to-scaling exponent 
Ai is still fairly large. It would therefore be of great interest to perform high-precision 
simulations at these values of g, going to very high values of L. 
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Table 2: Exact values of the static critical exponents, rounded to five decimal places. 
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L 


q = 1.25 


q = 1.50 


q = 1.75 


q = 2.00 


16 


0.5011732(38) 


0.5025161(59) 


0.5040267(83) 


0.5056788(68) 


32 


0.5005202(20) 


0.5011648(32) 


0.5019472(47) 


0.5028427(41) 


64 


0.5002301(10) 


0.5005395(17) 


0.5009352(26) 


0.5014281(24) 


128 


0.5001022(5) 


0.5002494(9) 


0.5004485(14) 


0.5007133(14) 


256 


0.5000444(3) 


0.5001146(5) 


0.5002159(8) 


0.5003559(8) 


512 


0.5000199(2) 


0.5000530(3) 


0.5001046(6) 


0.5001763(6) 


1024 


0.5000089(2) 


0.5000243(3) 


0.5000510(5) 


0.5000884(6) 



L 


q = 2.25 


q = 2.50 


q = 2.75 


q = 3.00 


16 


0.5075026(91) 


0.5094773(117) 


0.5115651(149) 


0.5138241(136) 


32 


0.5038896(56) 


0.5050706(76) 


0.5064185(101) 


0.5078851(97) 


64 


0.5020129(34) 


0.5027090(49) 


0.5035416(68) 


0.5045052(68) 


128 


0.5010481(21) 


0.5014567(31) 


0.5019755(45) 


0.5025847(47) 


256 


0.5005396(12) 


0.5007837(19) 


0.5010878(29) 


0.5014847(33) 


512 


0.5002817(10) 


0.5004183(17) 


0.5006029(26) 


0.5008522(32) 


1021 


0.5001110(10) 


0.5002238(17) 


0.5003359(28) 


0.5001888(35) 



L 


q = 3.25 


q = 3.50 


q = 3.75 


q = 4.00 


16 


0.5161911(167) 


0.5186811(202) 


0.5213781(243) 


0.5240486(230) 


32 


0.5095378(124) 


0.5113335(157) 


0.5132856(197) 


0.5154722(196) 


64 


0.5056256(92) 


0.5069004(122) 


0.5083371(162) 


0.5100231(169) 


128 


0.5033402(67) 


0.5042207(95) 


0.5052769(133) 


0.5065189(148) 


256 


0.5019717(49) 


0.5025922(74) 


0.5033613(111) 


0.5043157(132) 


512 


0.5011757(51) 


0.5015847(82) 


0.5021351(128) 


0.5028279(170) 


1024 


0.5006989(59) 


0.5009736(100) 


0.5013665(173) 


0.5018837(242) 



Table 3: Static data for (A/"), from the Monte Carlo simulations at the critical point of 
the 2- dimensional random cluster model, as a function of L and q, obtained by combining 
the data for all available k values for each q. The quoted error bar corresponds to one 
standard deviation. 
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L 


q = 1.25 


q = 1.50 


q = 1.75 


q = 2.00 


16 


0.5533 ±0.0004 


1.0380 ±0.0006 


1.4902 ±0.0008 


1.9285 ±0.0006 


32 


0.6180 ±0.0004 


1.1940 ±0.0006 


1.7753 ±0.0009 


2.3779 ±0.0008 


64 


0.6672 ±0.0004 


1.3284 ±0.0007 


2.0348 ±0.0010 


2.8226 ±0.0010 


128 


0.7058 ±0.0004 


1.4420 ±0.0007 


2.2767 ±0.0012 


3.2651 ±0.0012 


256 


0.7360 ±0.0004 


1.5366 ±0.0007 


2.4957 ±0.0013 


3.7062 ±0.0014 


512 


0.7594 ±0.0006 


1.6200 ±0.0011 


2.7016 ±0.0020 


4.1456 ±0.0023 


1024 


0.7768 ±0.0010 


1.6911 ±0.0018 


2.8922 ±0.0035 


4.5918 ±0.0043 




L 


q = 2.25 


q = 2.50 


q = 2.75 


q = 3.00 


16 


2.3666 ±0.0008 


2.8023 ±0.0011 


3.2511 ±0.0014 


3.7039 ±0.0013 


32 


3.0210 ±0.0011 


3.7149 ±0.0016 


4.4747 ±0.0022 


5.3043 ±0.0022 


64 


3.7224 ±0.0015 


4.7651 ± 0.0022 


5.9700 ±0.0033 


7.3834 ±0.0035 


128 


4.4734 ±0.0019 


5.9631 ±0.0031 


7.8056 ± 0.0049 


10.1121 ±0.0057 


256 


5.2731 ±0.0024 


7.3397 ±0.0042 


10.0718 ±0.0072 


13.6991 ±0.0092 


512 


6.1406 ±0.0043 


8.9200 ±0.0081 


12.8590 ±0.0141 


18.4085 ±0.0206 


1021 


7.0693 ± 0.0085 


10.7582 ±0.0171 


10.2779 ± 0.031-1 


2 i. 6809 ± 0.0521 



L 


q = 3.25 


q = 3.50 


q = 3.75 


q = 4.00 


16 


4.1724 ±0.0017 


4.6485 ± 0.0022 


5.1301 ±0.0027 


5.6282 ±0.0027 


32 


6.1993 ±0.0030 


7.1820 ±0.0040 


8.2551 ±0.0054 


9.4058 ±0.0058 


64 


9.0259 ±0.0052 


10.9483 ±0.0076 


13.1591 ±0.0111 


15.6961 ±0.0128 


128 


12.9732 ±0.0091 


16.5492 ±0.0144 


20.9891 ±0.0229 


26.4882 ±0.0294 


256 


18.5468 ±0.0158 


24.9710 ±0.0274 


33.5950 ±0.0481 


45.0709 ±0.0680 


512 


26.3990 ±0.0392 


37.7809 ±0.0751 


54.1813 ±0.1423 


78.3598 ±0.2298 


1024 


37.3057 ±0.1063 


56.9881 ± 0.2229 


87.9462 ± 0.4814 


136.9444 ±0.8705 



Table 4: Static data for Ch, from the Monte Carlo simulations at the critical point of the 
2-dimensional random cluster model, as a function of L and q, obtained by combining the 
data for all available k values for each pair {q,L). The quoted error bar corresponds to 
one standard deviation. 
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L 


q = 1.25 


g= 1.50 


q = 1.75 


q = 2.00 


16 


145.2585 ±0.0075 


142.6486±0.0106 


140.8411 ±0.0140 


139.6098 ±0.0107 


32 


497.8233 ±0.0265 


485.2393 ±0.0382 


476.2855 ±0.0515 


469.9356 ±0.0405 


64 


1706.4867 ±0.0932 


1650.6023 ±0.1357 


1610.1083 ±0.1860 


1581.4235 ±0.1495 


128 


5850.3427 ±0.3249 


5614.3344 ±0.4780 


5442.1190 ±0.6643 


5319.5479 ±0.5429 


256 


20054.6950 ±1.1280 


19099.3800 ±1.6687 


18399.0050 ±2.3542 


17890.2105 ±1.9519 


512 


68770.5920 ±5.3616 


64981.2580 ±7.9347 


62211.4880 ±11.2476 


60146.9825 ±9.8673 


1024 


235766.3800 ±30.2818 


221047.8000 ±45.1926 


210455.2400 ±64.7251 


202355.8500 ±56.1135 



L 


q = 2.25 


q = 2.50 


q = 2.75 


q = 3.00 


16 


138.8711 ±0.0134 


138.5052 ±0.0163 


138.3926 ±0.0195 


138.5876 ±0.0170 


32 


465.8207 ±0.0515 


463.3878 ±0.0642 


462.5502 ±0.0791 


462.6740 ±0.0708 


64 


1561.0853±0.1941 


1548.2598 ±0.2482 


1542.4292 ±0.3129 


1541.8545 ±0.2886 


128 


5232.5878 ±0.7205 


5173.4167 ±0.9438 


5144.4281 ±1.2191 


5134.3089 ±1.1567 


256 


17520.7785 ±2.6405 


17280.0615 ±3.5374 


17129.3300 ±4.7118 


17087.6517 ±4.5949 


512 


58713.4455 ±13.4485 


57681.1135 ±18.6703 


57059.7820 ±23.9577 


56845.9733 ±25.6428 


1024 


196609.0350 ±78.3919 


192554.3100 ±109.9204 


190128.7000 ±153.6076 


189071.0300 ±159.4872 



L 


q = 3.25 


q = 3.50 


q = 3.75 


q = 4.00 


16 


138.9566 ±0.0199 


139.5084 ±0.0231 


140.3091 ±0.0266 


141.0928 ±0.0242 


32 


464.1562 ±0.0851 


466.4774 ±0.1019 


469.6013 ±0.1210 


473.9444 ±0.1140 


64 


1546.5943 ±0.3570 


1555.8520 ±0.4412 


1569.0339 ±0.5429 


1588.4126 ±0.5280 


128 


5150.8298 ±1.4817 


5182.9838 ±1.8976 


5238.2491 ±2.4158 


5312.8763±2.4510 


256 


17123.9467 ±6.0699 


17256.6300 ±8.0307 


17479.6003 ±10.6954 


17798.6190 ±11.3951 


512 


56959.4603 ±35.1990 


57362. 2313±48.4641 


58221. 5237±65.2645 


59402.9003 ±75.5415 



1024 189355.8433 ±223.7014 190717.4067 ±322.6179 194007.2267 ±466.2996 198760.7775 ±550.8406 



Table 5: Static data for x, from the Monte Carlo simulations at the critical point of the 
2-dimensional random cluster model, as a function of L and q, obtained by combining the 
data for all available k values for each pair {q,L). The quoted error bar corresponds to 
one standard deviation. 
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L 


q = 1.25 


q = 1.50 


q = 1.75 


q = 2.00 


16 


0.92329(10) 


0.91489(13) 


0.91142(17) 


0.91137(13) 


32 


0.91903(10) 


0.91064(13) 


0.90705(17) 


0.90702(13) 


64 


0.91784(10) 


0.90927(14) 


0.90559(18) 


0.90596(14) 


128 


0.91755(10) 


0.90869(14) 


0.90492(19) 


0.90553(15) 


256 


0.91727(10) 


0.90868(14) 


0.90482(20) 


0.90511(16) 


512 


0.91748(14) 


0.90867(20) 


0.90500(27) 


0.90433(24) 


1024 


0.91729(23) 


0.90855(33) 


0.90595(46) 


0.90467(40) 




L 


q = 2.25 


q = 2.50 


q = 2.75 


q = 3.00 


16 


0.91427(16) 


0.91952(19) 


0.92604(22) 


0.93484(20) 


32 


0.91015(17) 


0.91564(21) 


0.92363(26) 


0.93262(23) 


64 


0.90856(18) 


0.91384(23) 


0.92198(29) 


0.93200(27) 


128 


0.90866(20) 


0.91375(26) 


0.92226(33) 


0.93190(32) 


256 


0.90774(21) 


0.91372(29) 


0.92091(38) 


0.93187(37) 


512 


0.90817(32) 


0.91313(45) 


0.92086(57) 


0.93186(61) 


1021 


0.90777(56) 


0.91280(78) 


0.92117(109) 


0.93191(114) 



L 


q = 3.25 


q = 3.50 


q = 3.75 


q = 4.00 


16 


0.94453(23) 


0.95554(27) 


0.96872(31) 


0.98129(29) 


32 


0.94431(28) 


0.95731(33) 


0.97163(40) 


0.98913(38) 


64 


0.94428(33) 


0.95874(41) 


0.97494(51) 


0.99544(50) 


128 


0.94522(40) 


0.95967(52) 


0.97789(66) 


0.99884(67) 


256 


0.94465(49) 


0.96087(65) 


0.98047(86) 


1.00431(92) 


512 


0.94516(84) 


0.96028(115) 


0.98054(155) 


1.00374(178) 


1024 


0.94561(160) 


0.95999(229) 


0.98103(329) 


1.00762(386) 



Table 6: Static data for ^/L, from the Monte Carlo simulations at the critical point of 
the 2-dimensional random cluster model, as a function of L and q, obtained by combining 
the data for all available k values for each q. The quoted error bar corresponds to one 
standard deviation. 
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L 


g= 1.25 


g= 1.50 


q = 1.75 


q = 2.00 


16 


186.2505 ±0.0063 


184.4023 ±0.0090 


183.1015±0.0119 


182.1899 ±0.0092 


32 


689.6841 ±0.0241 


680.2851 ±0.0349 


673.4366 ±0.0472 


668.4742 ±0.0373 


64 


2554.0565 ±0.0912 


2509.5209 ±0.1337 


2476.5244 ±0.1844 


2452.6604 ±0.1491 


128 


9458.3576 ±0.3426 


9256.6841 ±0.5088 


9106.1186 ±0.7136 


8996.7788 ±0.5878 


256 


35023.7950 ±1.2827 


34147.6420 ±1.9207 


33487.5700 ±2.7425 


32997.3910 ±2.2970 


512 


129717.7400 ±6.5774 


125972.7400 ±9.8879 


123159.5900 ±14.2165 


121001.7950 ±12.6325 


1024 


480361.5500 ±40.1000 


464676.1500 ±60.9749 


453086.4700 ±88.8330 


443892.0650 ±78.1385 



L 


q = 2.25 


q = 2.50 


q = 2.75 


q = 3.00 


16 


181.6159±0.0114 


181.3067 ±0.0140 


181.1556 ±0.0168 


181.2307 ±0.0147 


32 


665.1430 ±0.0477 


663.0250 ±0.0597 


662.0820 ±0.0738 


661.7617 ±0.0663 


64 


2435.0770 ±0.1948 


2423.3282 ±0.2505 


2417.2209 ±0.3176 


2415.0313±0.2944 


128 


8916.5343 ±0.7865 


8859.1099 ±1.0381 


8828.1239 ±1.3495 


8811.6847 ±1.2878 


256 


32630.0710 ±3.1391 


32381.5370 ±4.2429 


32211.1000 ±5.6996 


32145.2913 ±5.5950 


512 


119466.6650 ±17.4240 


118314.9800 ±24.4583 


117572.5400 ±31.6963 


117255.6470 ±34.1426 


1024 


437205.8550 ±110.7533 


432342.1300 ±157.0958 


429266.3250 ±222.0290 


427635.5230 ±232.5611 



L 


q = 3.25 


q = 3.50 


q = 3.75 


q = 4.00 


16 


181.4135 ±0.0172 


181.7222 ±0.0200 


182.2158 ±0.0231 


182.6714±0.0212 


32 


662.5142 ±0.0800 


663.8117 ±0.0962 


665.6519 ±0.1145 


668.4296 ±0.1084 


64 


2417.0586±0.3658 


2422.3096 ±0.4545 


2430.5493 ±0.5611 


2443.8329 ±0.5482 


128 


8818.7426 ±1.6581 


8836.6149 ±2.1356 


8874.0274 ±2.7312 


8926.2274 ±2.7887 


256 


32146.4423 ± 7.4382 


32235.6990 ±9.8938 


32400.5450 ± 13.2437 


32647.2718 ±14.1668 


512 


117238.6570 ±47.2276 


117484.8100 ±65.4957 


118192.9730 ±88.6774 


119114.2730 ±103.0354 


1024 


427561.8030 ±328.2653 


428396.6700 ±476.5507 


431327.5100±691.8068 


435425.1580 ±822.7561 



Table 7: Static data for Ci, from the Monte Carlo simulations at the critical point of the 
2-dimensional random cluster model, as a function of L and g, obtained by combining 
the data for all available k values for each q. The quoted error bar corresponds to one 
standard deviation. 
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L 


q = 1.25 


q = 1.50 


q = 1.75 


q = 2.00 


16 


1.3601 ±0.0008 


2.6552 ±0.0022 


4.5086 ±0.0048 


7.0555 ±0.0093 


32 


1.4516 ±0.0009 


2.9749 ±0.0026 


5.2959 ±0.0060 


8.6748 ±0.0127 


64 


1.5433 ±0.0010 


3.2988 ±0.0030 


6.1255 ±0.0075 


10.4931 ±0.0167 


128 


1.6225 ±0.0011 


3.6209 ±0.0034 


7.0060 ±0.0093 


12.4925 ±0.0217 


256 


1.7014±0.0012 


3.9284 ±0.0039 


7.8884 ±0.0110 


14.7173 ±0.0278 


512 


1.7679 ±0.0016 


4.2341 ±0.0060 


8.8097 ±0.0176 


17.0928 ±0.0493 


1024 


1.8281 ±0.0028 


4.5295 ±0.0109 


9.8139 ±0.0340 


19.8144 ±0.0982 



L 


q = 2.25 


q = 2.50 


q = 2.75 


q = 3.00 


16 


10.5266 ±0.0169 


15.0260 ±0.0288 


20.8370 ±0.0469 


28.2660 ± 0.0739 


32 


13.4708 ±0.0243 


20.0908 ±0.0443 


29.2540 ±0.0778 


41.3643 ±0.1309 


64 


17.0098 ±0.0346 


26.6206 ±0.0675 


40.4635 ±0.1264 


60.2637 ±0.2298 


128 


21.2642 ±0.0483 


34.7971 ±0.1008 


55.3810 ±0.2025 


87.1456 ±0.3990 


256 


26.1625 ±0.0657 


44.9461 ±0.1479 


75.4692 ±0.3214 


124.4389 ±0.6808 


512 


32.0642 ±0.1234 


57.4369 ±0.3026 


101.4174 ±0.7096 


175.6975 ±1.6194 


1024 


38.0893 ±().2()o7 


73.2071 ± 0.0912 


13o.0123± 1.7390 


240.8774 ±4,2833 



L 


q = 3.25 


q = 3.50 


q = 3.75 


q = 4.00 


16 


37.4020 ±0.1125 


48.7655 ±0.1673 


62.7027 ±0.2440 


79.4340 ±0.3475 


32 


57.1966 ±0.2126 


78.7249 ±0.3431 


106.2816 ±0.5379 


142.0782 ±0.8310 


64 


87.3178 ±0.4006 


126.1715 ±0.6960 


181.3644 ±1.1991 


253.4333 ±1.9791 


128 


134.4255 ±0.7644 


206.2906 ±1.4542 


312.1399 ±2.7062 


462.6885 ±4.8818 


256 


202.0761 ±1.4086 


329.0069 ±2.9261 


534.9698 ±6.0685 


875.8576 ±12.7193 


512 


309.5691 ±3.7861 


532.7875 ±8.5469 


933.1742 ±19.8120 


1652.8786 ±46.7469 


1024 


446.2539 ± 10.4070 


822.7567 ±26.0508 


1581.9426 ±69.4509 


3049.0566 ±186.3115 



Table 8: Estimates of Tint,5' from the Monte Carlo simulations at the critical point of the 
2-dimensional random cluster model, as a function of L and q, ior k — 1. The quoted 
error bar corresponds to one standard deviation. 
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L 


q = 1.25 


q = 1.50 


q = 1.75 


q = 2.00 


16 


1.2040 ±0.0007 


2.2847 ±0.0019 


3.8600 ±0.0041 


6.0689 ±0.0080 


32 


1.2656 ±0.0008 


2.5269 ±0.0022 


4.5125 ±0.0052 


7.4842 ±0.0110 


64 


1.3202 ±0.0009 


2.7598 ±0.0025 


5.1783 ±0.0063 


9.0496 ±0.0144 


128 


1.3614 ±0.0009 


2.9810 ±0.0028 


5.8732 ±0.0078 


10.7689 ±0.0187 


256 


1.4015 ±0.0010 


3.1840 ±0.0032 


6.5651 ±0.0092 


12.6943 ±0.0240 


512 


1.4325 ±0.0013 


3.3804 ±0.0048 


7.2791 ±0.0146 


14.7495 ±0.0425 


1024 


1.4560 ±0.0022 


3.5685 ±0.0086 


8.0614 ±0.0279 


17.1123 ±0.0848 



L 


q = 2.25 


q = 2.50 


q = 2.75 


q = 3.00 


16 


9.1341 ±0.0147 


13.1592 ±0.0252 


18.4381 ±0.0415 


25.2445 ±0.0660 


32 


11.7916±0.0213 


17.8549 ±0.0394 


26.3710 ±0.0701 


37.7574±0.1195 


64 


15.0067 ±0.0306 


23.9303 ±0.0607 


37.0129±0.1156 


55.9353 ±0.2133 


128 


18.8651 ±0.0428 


31.6032 ±0.0916 


51.3588 ±0.1878 


82.0584 ±0.3757 


256 


23.3330 ±0.0586 


41.2200 ±0.1356 


70.7459 ±0.3013 


118.5736 ±0.6487 


512 


28.7842±0.1107 


53.1546 ±0.2800 


96.0301 ±0.6719 


169.0612 ±1.5582 


1024 


34.9401 ±0.2400 


G8.3955 ± ()Mo2 


129.5483 ±l.(i(il8 


239.4339 ± 4. l.')41 



L 


q = 3.25 


q = 3.50 


q = 3.75 


q = 4.00 


16 


33.6892±0.1013 


44.2837±0.1520 


57.2882 ±0.2229 


73.0397±0.3195 


32 


52.7843 ±0.1962 


73.3099 ±0.3195 


99.7789 ±0.5050 


134.2543 ±0.7852 


64 


82.0705 ±0.3765 


119.8115 ±0.6609 


173.5483 ±1.1474 


244.1495 ±1.9066 


128 


128.1915±0.7289 


198.6992 ±1.4006 


302.8677±2.6258 


451.5682 ±4.7645 


256 


194.9751 ±1.3591 


320.3198 ±2.8488 


524.1844 ±5.9461 


862.6987 ±12.5282 


512 


301.2867 ±3.6848 


522.8481 ±8.3875 


921.2799 ±19.5595 


1638.0215 ±46.3267 


1024 


437.3916 ±10.2003 


812.1417 ±25.7147 


1568.8474 ±68.8760 


3032.3859 ±185.2928 



Table 9: Estimates of Tint,A/' from the Monte Carlo simulations at the critical point of the 
2-dimensional random cluster model, as a function of L and q, ior k — 1. The quoted 
error bar corresponds to one standard deviation. 
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L 


q = 1.25 


g= 1.50 


g= 1.75 


q = 2.00 


16 


1.3386 ±0.0008 


2.6034 ±0.0021 


4.4205 ±0.0047 


6.9286 ±0.0091 


32 


1.4078 ±0.0009 


2.8602 ±0.0025 


5.0735 ±0.0058 


8.3286 ±0.0122 


64 


1.4676 ±0.0010 


3.0892 ±0.0028 


5.7122 ±0.0070 


9.8003 ±0.0156 


128 


1.5118±0.0010 


3.2941 ±0.0031 


6.3331 ±0.0084 


11.2921 ±0.0196 


256 


1.5495±0.0011 


3.4664 ±0.0034 


6.8843 ±0.0096 


12.8477 ±0.0243 


512 


1.5770 ±0.0015 


3.6264 ±0.0051 


7.4316 ±0.0149 


14.3819 ±0.0415 


1024 


1.6009 ±0.0024 


3.7573 ±0.0090 


7.9652 ±0.0276 


16.0737 ±0.0796 



L 


q = 2.25 


q = 2.50 


q = 2.75 


q = 3.00 


16 


10.3488 ±0.0166 


14.8081 ±0.0284 


20.5896 ±0.0464 


28.0138 ±0.0733 


32 


12.9917±0.0235 


19.4472 ±0.0429 


28.4240 ±0.0756 


40.4259 ±0.1279 


64 


15.9573 ±0.0325 


25.1767 ±0.0639 


38.5615 ±0.1205 


57.9953 ±0.2212 


128 


19.3490 ±0.0439 


32.1053 ±0.0930 


51.5394 ±0.1884 


82.2076 ±0.3764 


256 


23.0132 ±0.0578 


40.1938 ±0.1322 


68.5937±0.2922 


115.2704 ±0.6306 


512 


27.2615 ±0.1049 


49.7110 ±0.2619 


89.9306 ±0.6292 


159.7078 ±1.4720 


1024 


31.7804 ±0.2183 


01.4.579 ±0.3798 


110.7933 ±1.4982 


220.0990 ±3.8291 



L 


q = 3.25 


q = 3.50 


q = 3.75 


q = 4.00 


16 


37.2039±0.1119 


48.6233 ±0.1669 


62.6624 ±0.2439 


79.6370 ±0.3484 


32 


56.2252 ±0.2090 


77.8100 ±0.3391 


105.4646 ±0.5338 


141.5583 ±0.8279 


64 


84.5549 ±0.3879 


123.3440 ±0.6804 


178.5750 ±1.1806 


250.9246 ±1.9595 


128 


128.5492 ±0.7310 


199.5134 ±1.4064 


305.0319 ±2.6446 


457.0131 ±4.8219 


256 


190.1111 ±1.3252 


315.2552 ±2.8038 


520.8366 ±5.9082 


860.9974 ±12.5035 


512 


287.5530 ±3.5168 


506.8770 ±8.1313 


898.2129 ±19.0697 


1617.6331 ±45.7501 


1024 


408.0224 ±9.5154 


777.6178 ±24.6216 


1526.7086 ±67.0260 


2980.4947 ±182.1220 



Table 10: Estimates of rint,52 from the Monte Carlo simulations at the critical point of 
the 2-dimensional random cluster model, as a function of L and q, ior k — 1. The quoted 
error bar corresponds to one standard deviation. 
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L 


g= 1.25 


g= 1.50 


g= 1.75 


q = 2.00 


16 


1.1821 ±0.0007 


2.1678 ±0.0018 


3.5450 ±0.0038 


5.4257 ±0.0072 


32 


1.2066 ±0.0007 


2.2716 ±0.0020 


3.8255 ±0.0044 


6.0790 ±0.0089 


64 


1.2285 ±0.0008 


2.3601 ±0.0021 


4.1069 ±0.0050 


6.7601 ±0.0108 


128 


1.2433 ±0.0008 


2.4422 ±0.0023 


4.3807 ±0.0058 


7.4396 ±0.0129 


256 


1.2563 ±0.0009 


2.5088 ±0.0025 


4.6063 ±0.0064 


8.1526±0.0154 


512 


1.2672 ±0.0012 


2.5784 ±0.0036 


4.8445 ±0.0097 


8.8373 ±0.0255 


1024 


1.2756 ±0.0019 


2.6249 ±0.0063 


5.0706 ±0.0176 


9.5914 ±0.0475 



L 


q = 2.25 


q = 2.50 


q = 2.75 


q = 3.00 


16 


7.9521 ±0.0128 


11.2607 ±0.0216 


15.5512 ±0.0350 


21.1189 ±0.0552 


32 


9.2758 ±0.0168 


13.7011 ±0.0302 


19.8378 ±0.0528 


28.2005 ±0.0892 


64 


10.7242 ±0.0218 


16.6498 ±0.0422 


25.3352 ±0.0791 


38.0837 ±0.1452 


128 


12.3646 ±0.0281 


20.2411 ±0.0587 


32.2903 ±0.1181 


51.5226 ±0.2359 


256 


14.1174 ±0.0354 


24.2638 ±0.0798 


41.3352 ±0.1761 


70.0114 ±0.3830 


512 


16.1705 ±0.0622 


29.0610 ±0.1531 


52.5554 ±0.3677 


94.2286 ±0.8685 


1024 


18. .380() ±0.1202 


34.7041 ±0.3274 


GO. 1941 ±0.8491 


128.1428 ±2.2233 



L 


q = 3.25 


q = 3.50 


q = 3.75 


q = 4.00 


16 


28.0987 ±0.0845 


36.7802 ±0.1262 


47.5330 ±0.1850 


60.6406 ±0.2653 


32 


39.3515 ±0.1463 


54.7327 ±0.2386 


74.6290 ±0.3777 


100.4741 ±0.5877 


64 


55.6487 ±0.2553 


82.1045 ±0.4529 


119.8880 ±0.7926 


170.2119 ±1.3292 


128 


81.4775 ±0.4633 


127.4713 ±0.8986 


197.8031 ±1.7149 


300.7931 ±3.1737 


256 


116.2497 ±0.8103 


196.8057 ±1.7503 


330.7032 ±3.7514 


557.2872 ±8.0930 


512 


171.9805 ±2.1034 


308.7344 ±4.9527 


557.6357 ±11.8390 


1014.8762 ±28.7029 


1024 


239.5051 ±5.5854 


471.3493 ±14.9243 


943.8663 ±41.4379 


1849.8784 ±113.0361 



Table 11: Estimates of rint,^:-' from the Monte Carlo simulations at the critical point of 
the 2-dimensional random cluster model, as a function of L and q, for k — 1. The quoted 
error bar corresponds to one standard deviation. 
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L 


q = 1.25 


g= 1.50 


g= 1.75 


q = 2.00 


16 


1.3101 ±0.0008 


2.5278 ±0.0021 


4.2775 ±0.0046 


6.6989 ±0.0088 


32 


1.3708 ±0.0008 


2.7571 ±0.0024 


4.8684 ±0.0056 


7.9857±0.0117 


64 


1.4224 ±0.0009 


2.9585 ±0.0027 


5.4412 ±0.0067 


9.3266 ±0.0149 


128 


1.4596 ±0.0009 


3.1373 ±0.0030 


5.9956 ±0.0079 


10.6786 ±0.0186 


256 


1.4918±0.0010 


3.2868 ±0.0033 


6.4856 ±0.0091 


12.0813 ±0.0229 


512 


1.5155 ±0.0014 


3.4287 ±0.0048 


6.9711 ±0.0139 


13.4617 ±0.0388 


1024 


1.5364 ±0.0023 


3.5411 ±0.0085 


7.4506 ±0.0258 


14.9718 ±0.0742 



L 


q = 2.25 


q = 2.50 


q = 2.75 


q = 3.00 


16 


10.0086 ±0.0161 


14.3476 ±0.0275 


20.0130 ±0.0451 


27.2880 ±0.0714 


32 


12.4722 ±0.0225 


18.7396 ±0.0414 


27.4966 ±0.0731 


39.3047±0.1244 


64 


15.2204 ±0.0310 


24.0900 ±0.0611 


37.1348 ±0.1160 


56.1933 ±0.2143 


128 


18.3474 ±0.0416 


30.5895 ±0.0886 


49.4474 ±0.1808 


79.3969 ±0.3635 


256 


21.6896 ±0.0544 


38.1134±0.1254 


65.5494 ±0.2792 


111.0562 ±0.6076 


512 


25.5866 ±0.0984 


46.9593 ±0.2474 


85.7120 ±0.5997 


153.1954 ±1.4120 


1024 


29.7135 ±0.2041 


57.6989 ±0.3443 


110.7213 ±1.4203 


211.5007 ± 3. (i()95 



L 


q = 3.25 


q = 3.50 


q = 3.75 


q = 4.00 


16 


36.4019 ±0.1095 


47.7804 ±0.1640 


61.9085 ±0.2409 


79.1044 ±0.3461 


32 


54.9946 ±0.2045 


76.5469 ±0.3336 


104.2901 ±0.5278 


141.5283 ±0.8278 


64 


82.5872 ±0.3789 


121.5950 ±0.6708 


177.0526 ±1.1706 


251.4255 ±1.9634 


128 


125.2155±0.7120 


196.1776 ±1.3829 


303.5483 ±2.6317 


461.4449±4.8687 


256 


185.1345 ±1.2905 


310.3888 ±2.7605 


518.9183 ±5.8864 


868.1885 ±12.6080 


512 


280.0361 ±3.4249 


500.1480 ±8.0233 


900.8537 ±19.1258 


1637.5816 ±46.3143 


1024 


395.2831 ±9.2183 


762.0676 ±24.1293 


1519.4102 ±66.7056 


3054.3411 ±186.6344 



Table 12: Estimates of Tint,Ci from the Monte Carlo simulations at the critical point of 
the 2-dimensional random cluster model, as a function of L and q, ior k — 1. The quoted 
error bar corresponds to one standard deviation. 
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L 


q = 1.25 


q = 1.50 


q = 1.75 


q = 2.00 


16 


1.0547 ±0.0007 


1.8490 ±0.0015 


2.9362 ±0.0031 


4.3973 ±0.0058 


32 


1.0744 ±0.0007 


1.9293±0.0017 


3.1495 ±0.0036 


4.8826 ±0.0071 


64 


1.0920 ±0.0007 


1.9968 ±0.0018 


3.3571 ±0.0041 


5.3808 ±0.0086 


128 


1.1020 ±0.0007 


2.0563 ±0.0020 


3.5553 ±0.0047 


5.8582 ±0.0102 


256 


1.1120±0.0008 


2.1031 ±0.0021 


3.7110 ±0.0052 


6.3568 ±0.0120 


512 


1.1195 ±0.0010 


2.1560 ±0.0030 


3.8788 ±0.0078 


6.8290 ±0.0197 


1024 


1.1265±0.0017 


2.1835 ±0.0052 


4.0403 ±0.0140 


7.3298 ±0.0363 



L 


q = 2.25 


g = 2.50 


q = 2.75 


q = 3.00 


16 


6.3209 ±0.0102 


8.8061 ±0.0169 


12.0317 ±0.0271 


16.1764 ±0.0423 


32 


7.3082 ±0.0132 


10.6241 ±0.0234 


15.2020 ±0.0404 


21.3897±0.0677 


64 


8.3609 ±0.0170 


12.7509 ±0.0323 


19.1601 ±0.0599 


28.6341 ±0.1092 


128 


9.5115 ±0.0216 


15.2680 ±0.0442 


24.1050 ±0.0881 


38.1520 ±0.1747 


256 


10.7245 ±0.0269 


18.0907 ±0.0595 


30.4051 ±0.1295 


51.3142 ±0.2807 


512 


12.1338 ±0.0467 


21.4336 ±0.1129 


38.3258 ±0.2682 


68.1149 ±0.6278 


1024 


13.0414 ±0.0937 


23.2171 ±0.2379 


47.5954 ±0.(jlU5 


92.1484 ± 1.3988 



L 


q = 3.25 


q = 3.50 


q = 3.75 


q = 4.00 


16 


21.3907±0.0643 


27.8259 ±0.0955 


35.8844±0.1397 


45.7629 ±0.2002 


32 


29.7129 ±0.1105 


41.2669 ±0.1799 


56.1696 ±0.2843 


76.2426 ±0.4459 


64 


41.7815 ±0.1917 


61.7907 ±0.3409 


90.4411 ±0.5979 


129.5723 ±1.0118 


128 


60.2645 ±0.3427 


95.1785 ±0.6709 


148. 7277± 1.2894 


231.9678 ±2.4475 


256 


85.4706 ±0.5958 


145.8670 ±1.2973 


249.9941 ±2.8358 


429.9691 ±6.2441 


512 


125.2840 ±1.5322 


229.6477 ±3.6840 


425.3373 ±9.0302 


810.2299 ±22.9150 


1024 


174.4576 ±4.0685 


348.2103 ±11.0253 


724.6440 ±31.8135 


1487.0973 ±90.8685 



Table 13: Estimates of rint,C2 from the Monte Carlo simulations at the critical point of 
the 2-dimensional random cluster model, as a function of L and q, ior k — 1. The quoted 
error bar corresponds to one standard deviation. 
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L 


g= 1.25 


g= 1.50 


g= 1.75 


q = 2.00 


16 


1.0627 ±0.0007 


1.9037 ±0.0015 


3.0988 ±0.0033 


4.7417 ±0.0063 


32 


1.0882 ±0.0007 


2.0070 ±0.0017 


3.3763 ±0.0039 


5.3802 ±0.0079 


64 


1.1096 ±0.0007 


2.0948 ±0.0019 


3.6466 ±0.0045 


6.0321 ±0.0096 


128 


1.1247 ±0.0007 


2.1761 ±0.0021 


3.9002 ±0.0051 


6.6878 ±0.0116 


256 


1.1379 ±0.0008 


2.2400 ±0.0022 


4.1293 ±0.0058 


7.3441 ±0.0139 


512 


1.1487 ±0.0011 


2.3012 ±0.0032 


4.3459 ±0.0087 


7.9879 ±0.0230 


1024 


1.1573 ±0.0018 


2.3518 ±0.0056 


4.5466 ±0.0158 


8.6848 ±0.0430 



L 


q = 2.25 


q = 2.50 


q = 2.75 


q = 3.00 


16 


6.9598 ±0.0112 


9.8960 ±0.0190 


13.7462 ±0.0310 


18.6750 ±0.0488 


32 


8.2349 ±0.0149 


12.2476 ±0.0270 


17.8387 ±0.0474 


25.5572 ±0.0809 


64 


9.6298 ±0.0196 


15.0431 ±0.0382 


23.1057 ±0.0722 


35.0341 ±0.1336 


128 


11.1923 ±0.0254 


18.4369 ±0.0534 


29.6671 ±0.1085 


47.6201 ±0.2180 


256 


12.8179 ±0.0322 


22.2087 ±0.0731 


38.1085 ±0.1623 


65.1028 ±0.3562 


512 


14.73.38 ±0.0567 


26.6674 ±0.1405 


48.6601 ±0.3405 


87.6815 ±0.8081 


1024 


1G.7174± 0.1148 


31.<sr49± 0.3007 


01.3270 ±0.7892 


118.9944 ± 2. 0()43 



L 


q = 3.25 


q = 3.50 


q = 3.75 


q = 4.00 


16 


24.9741 ±0.0751 


32.8814±0.1128 


42.7608 ±0.1664 


54.9178 ±0.2403 


32 


35.7876 ±0.1330 


50.1792 ±0.2187 


68.8743 ±0.3486 


94.4193 ±0.5522 


64 


51.7567 ±0.2374 


77.0195 ±0.4249 


113.3164 ±0.7492 


163.7752 ±1.2789 


128 


76.0751 ±0.4326 


120. 8170±0. 8516 


190.4901 ±1.6515 


297.1283 ±3.1350 


256 


109.5690 ±0.7638 


188.0194 ±1.6722 


321.0373 ±3.6417 


553.7972 ±8.0423 


512 


162.8091 ±1.9912 


298.8569 ±4.7942 


551.3396 ±11.7054 


1044.3142 ±29.5354 


1024 


226.4912 ±5.2819 


452.2253 ±14.3188 


928.5042 ±40.7635 


1950.0423 ±119.1566 



Table 14: Estimates of Tint.Cs from the Monte Carlo simulations at the critical point of 
the 2-dimensional random cluster model, as a function of L and q, ior k — 1. The quoted 
error bar corresponds to one standard deviation. 
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Fits to ilL = x' 








Fits to C/L = 


X* + BL-P 










q 


-^min 


X* 




DF CL(%) 


^min 


X* 


B 


P 




DF CL(%) 


Ai (exact) 


1.25 


256 


0.91734(8) 


1.48 


2 


48 


16 


0.91735(7) 


0.900(243) 


1.81(10) 


2.53 


4 


64 


1.80702 


1.50 


128 


0.90867(9) 


0.17 


3 


98 


16 


0.90857(10) 


0.563(173) 


1.62(11) 


0.77 


4 


94 


1.63551 


1.75 


128 


0.90496(12) 


5.19 


3 


16 


16 


0.90487(14) 


0.582(228) 


1.62(14) 


6.93 


4 


14 


1.47905 


2.00 


512 


0.90442(20) 


0.54 


1 


46 


32 


0.90387(74) 


0.017(9) 


0.49(21) 


4.12 


3 


25 


1.33333 


2.25 


256 


0.90786(17) 


1.30 


2 


52 


[16] 


0.90799(16) 


0.465(184) 


1.55(15) 


8.19 


4 


8 


1.19517 


2.50 


64 


0.91369(14) 


3.36 


4 


50 


16 


0.91335(21) 


0.396(195) 


1.50(18) 


3.69 


4 


45 


1.06191 


2.75 


256 


0.92092(30) 


0.06 


2 


97 


16 


0.92085(45) 


0.066(33) 


0.92(20) 


5.96 


4 


20 


0.93110 


3.00 


64 


0.93193(17) 


0.13 


4 


100 


16 


0.93184(23) 


0.704(1.120) 1.97(59) 


0.02 


4 


100 


0.80000 


3.25 


128 


0.94503(29) 


0.97 


3 


81 
















0.66499 


3.50 


256 


0.96069(55) 


0.30 


2 


86 


16 


0.96223(198) 


-0.025(9) 


0.47(23) 


1.21 


4 


88 


0.51991 


3.75 


256 


0.98051(73) 


0.03 


2 


99 


32 


0.98591(540) 


-0.060(27) 


0.41(23) 


1.19 


3 


75 


0.34986 


4.00 


256 


1.00434(80) 


0.84 


2 


66 


16 


1.01323(335) 


-0.100(9) 


0.41(7) 


4.81 


4 


31 


0.00000 



Table 15: Fits for ^/L. The quoted error bar corresponds to one standard deviation. Lmin in brackets indicates a poor fit 
(confidence level < 10%). 





Fits to X = AL'^I'' 




q 






A 




DF CL(%) 


^/v (exact) 


deviation 


dev(o") 


1.25 


256 


1.77773(9) 


1.0496(6) 


1.40 


1 


24 


1.77764 


0.00009 


1.0 


1.50 


128 


1.76639(8) 


1.0645(5) 


0.30 


2 


86 


1.76644 


-0.00005 


-0.6 


1.75 


128 


1.75759(12) 


1.0768(7) 


2.34 


2 


31 


1.75738 


0.00022 


1.8 


2.00 


128 


1.74969(11) 


1.0937(6) 


2.54 


2 


28 


1.75000 


-0.00031 


-2.9 


2.25 


256 


1.74428(26) 


1.1039(17) 


1.37 


1 


24 


1.74404 


0.00024 


0.9 


2.50 


256 


1.73902(37) 


1.1209(24) 


0.01 


1 


94 


1.73932 


-0.00030 


-0.8 


2.75 


128 


1.73575(28) 


1.1316(17) 


0.85 


2 


65 


1.73576 


-0.00001 


-0.03 


3.00 


128 


1.73438(28) 


1.1371(17) 


0.72 


2 


70 


1.73333 


0.00105 


3.8 


3.25 


128 


1.73341(37) 


1.1460(23) 


0.30 


2 


86 


1.73214 


0.00127 


3.4 


3.50 


256 


1.73304(99) 


1.1571(67) 


0.01 


1 


93 


1.73246 


0.00058 


0.6 


3.75 


128 


1.73745(67) 


1.1431(40) 


1.46 


2 


48 


1.73517 


0.00229 


3.4 


4.00 


256 


1.73968(152) 


1.1502(101) 


0.44 


1 


51 


1.75000 


-0.01032 


-6.8 



Table 16: Fits of the susceptibility x to a pure power law. The quoted error bar corre- 
sponds to one standard deviation. Lmin in brackets indicates a poor fit (confidence level 
< 10%). The final two columns show the deviation of the estimated ^/v from the known 
exact value, in absolute terms and in units of its standard deviation. 
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Fits to x/L 




A 






Fits to x/L^^"" 


= A + BL- 


-p 












A 




DF CL(%) 


^min 


A 


B 


P 




DF CL(%) 


Ai (exact) 


1.25 


256 


1.05017(5) 


2.38 


2 


30 


16 


1.05013(8) 


0.011(5) 


0.87(19) 


4.72 


4 


32 


1.80702 


1.50 


128 


1.06417(6) 


0.70 


3 


87 


16 


1.06379(53) 


0.003(1) 


0.34(30) 


1.55 


4 


82 


1.63551 


1.75 


64 


1.07811(7) 


7.16 


4 


13 
















1.47905 


2.00 


512 


1.09148(15) 


0.50 


1 


48 


[16] 


1.09195(8) 


-0.231(316) 1.88(50) 


15.63 


4 


0.4 


1.33333 


2.25 


256 


1.10545(13) 


2.19 


2 


33 


[16] 


1.10565(15) 


-0.078(43) 


1.22(21) 


7.84 


4 


10 


1.19517 


2.50 


128 


1.11875(14) 


2.34 


3 


51 


16 


1.11912(26) 


-0.062(22) 


0.94(14) 


3.96 


4 


41 


1.06191 


2.75 


128 


1.13159(18) 


0.85 


3 


84 


16 


1.13192(29) 


-0.161(51) 


1.12(12) 


5.80 


4 


21 


0.93110 


3.00 


256 


1.14415(25) 


1.71 


2 


43 


16 


1.14541(49) 


-0.086(12) 


0.73(6) 


0.43 


4 


98 


0.80000 


3.25 


512 


1.15559(63) 


0.25 


1 


62 


16 


1.15656(67) 


-0.111(14) 


0.70(6) 


2.74 


4 


60 


0.66499 


3.50 


256 


1.16102(46) 


0.35 


2 


84 


16 


1.16365(96) 


-0.119(15) 


0.65(6) 


0.99 


4 


91 


0.51991 


3.75 


256 


1.15852(61) 


0.42 


2 


81 


16 


1.16358(199) 


-0.080(9) 


0.48(7) 


0.79 


4 


94 


0.34986 


4.00 


[512] 


1.07693(124) 


2.97 


1 


8 
















0.00000 



Table 17: Fits for x/L'^/'^ with 7/z/ set to its exact value fll.lOp . The quoted error bar corresponds to one standard deviation. 
Lmin in brackets indicates a poor fit (confidence level < 10%). 





Fits to Ci = AL'^-I^f" 




q 


^min 


/3/u 


A 




DF CL(%) 


P/v (exact) 


deviation 


dev((T) 


1.25 


256 


0.11111(6) 


0.9896(3) 


1.29 


1 


26 


0.11118 


-0.00007 


-1.2 


1.50 


16 


0.11681(2) 


0.9958(1) 


2.01 


5 


85 


0.11678 


0.00003 


1.6 


1.75 


128 


0.12115(8) 


1.0004(4) 


2.21 


2 


33 


0.12131 


-0.00016 


-2.1 


2.00 


128 


0.12518(7) 


1.0080(4) 


2.04 


2 


36 


0.12500 


0.00018 


2.7 


2.25 


256 


0.12787(17) 


1.0118(10) 


1.20 


1 


27 


0.12798 


-0.00011 


-0.6 


2.50 


128 


0.13033(13) 


1.0177(7) 


1.48 


2 


48 


0.13034 


-0.00001 


-0.10 


2.75 


128 


0.13230(18) 


1.0237(10) 


1.65 


2 


44 


0.13212 


0.00018 


1.0 


3.00 


64 


0.13283(11) 


1.0245(5) 


1.24 


3 


74 


0.13333 


-0.00050 


-4.6 


3.25 


128 


0.13366(24) 


1.0295(13) 


0.86 


2 


65 


0.13393 


-0.00027 


-1.1 


3.50 


64 


0.13315(19) 


1.0289(9) 


2.39 


3 


50 


0.13377 


-0.00062 


-3.3 


3.75 


32 


0.13174(15) 


1.0263(7) 


1.74 


4 


78 


0.13242 


-0.00067 


-4.4 


4.00 


256 


0.13200(104) 


1.0357(62) 


0.54 


1 


46 


0.12500 


0.00700 


6.8 



Table 18: Fits of the mean size Ci of the largest cluster to a pure power law. The quoted 

error bar corresponds to one standard deviation. Lmin in brackets indicates a poor fit 
(confidence level < 10%). The final two columns show the deviation of the estimated P/v 
from the known exact value, in absolute terms and in units of its standard deviation. 
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Fits to Ch = AL 


O-jv 






Fits to Ch = AL""/" + B 




q 




a/v 


A 




DF 


CL(%) 




a/v 


A 


B 




DF CL(%) 


a/v (exact) 


dev 


dcv((T) 


1.25 


[256] 


0.0406(9) 


0.588(3) 


16.40 


1 


5 X 10"^' 


16 


-0.3675(49) 


-0.791(6) 


0.839(2) 


1.82 


4 


77 


-0.35527 


-0.012 


-2.5 


1.50 


[256] 


0.0712(8) 


1.036(5) 


30.87 


1 


3 X 10^6 


128 


-0.1977(186) 


-1.933(23) 


2.182(59) 


0.37 


1 


54 


-0.22663 


0.029 


1.6 


1.75 


[256] 


0.1089(8) 


1.366(7) 


31.62 


1 


2 X 10-6 


128 


-0.0965(135) 


-5.413(347) 


5.665(440) 


0.15 


1 


70 


-0.10929 


0.013 


0.9 


2.00 


[256] 


0.1569(6) 


1.554(6) 


43.28 


1 


5 X 10-9 


(see text) 














O(log) 






2.25 


[256] 


0.2143(8) 


1.608(8) 


36.65 


1 


1 X lO-'^ 


128 


0.1108(65) 


5.856(558) 


-5.554(640) 


0.45 


1 


50 


0.10363 


0.007 


1.1 


2.50 


[256] 


0.2778(10) 


1.573(10) 


9.47 


1 


0.2 


64 


0.2025(32) 


3.421(109) 


-3.176(150) 


0.73 


2 


69 


0.20357 


-0.001 


-0.3 


2.75 


[256] 


0.3488(13) 


1.457(11) 


7.44 


1 


0.6 


32 


0.2995(18) 


2.293(33) 


-2.001(54) 


0.39 


3 


94 


0.30168 


-0.002 


-1.2 


3.00 


256 


0.4254(13) 


1.295(10) 


0.42 


1 


52 


32 


0.3947(15) 


1.681(18) 


-1.296(37) 


2.14 


3 


54 


0.40000 


-0.005 


-3.6 


3.25 


[256] 


0.5063(17) 


1.120(11) 


2.94 


1 


9 


32 


0.4902(16) 


1.275(14) 


-0.771(38) 


3.04 


3 


39 


0.50126 


-0.011 


-6.8 


3.50 


64 


0.5952(7) 


0.921(3) 


0.94 


3 


82 


64 


0.5945(35) 


0.925(22) 


-0.020(102) 


0.90 


2 


64 


0.61007 


-0.016 


-4.5 


3.75 


256 


0.6920(31) 


0.724(13) 


0.69 


1 


41 


128 


0.7095(86) 


0.635(38) 


1.149(369) 


0.03 


1 


86 


0.73760 


-0.028 


-3.3 


4.00 


256 


0.7998(35) 


0.534(11) 


0.35 


1 


55 


128 


0.8338(85) 


0.416(24) 


2.713(392) 


0.53 


1 


47 


1 (log-^/^) 


-0.166 


-19.6 



Table 19: Fits for the specific heat Ch- The quoted error bar corresponds to one standard deviation. Lmin in brackets indicates 
a poor fit (confidence level < 10%). The final two columns show the deviation of the estimated a/v (from the fit AL"^^^ + B) 
from the known exact value, in absolute terms and in units of its standard deviation. 



L 


q^2 


q = 3.25 


16 


2.143(3) 


2.020(7) 


32 


2.123(4) 


2.013(9) 


64 


2.100(4) 


1.990(11) 


128 


2.080(4) 


2.015(14) 


256 


2.069(5) 


1.991(17) 


512 


2.047(7) 


2.034(30) 


1024 


2.051(12) 


1.953(56) 



Table 20: Ratios Tint,£'{k — l)/Tint,e'{k = 2) as a function of L for g = 2 and q — 3.25. 
Error bars are one standard deviation. 







Fits to 


rmt,£'i^)/rint,s'{k) 


= A 




Q 


k 


-^min 


A 




DF 


CL(%) 


2.00 


2 


512 


2.048(6) 


0.09 


1 


76 


2.25 


2 


256 


2.036(5) 


0.51 


2 


77 


2.50 


2 


128 


2.022(5) 


0.02 


3 


100 


2.75 


2 


128 


2.006(6) 


0.87 


3 


83 


3.00 


2 


32 


2.015(5) 


3.63 


5 


60 


3.25 


2 


32 


2.005(6) 


5.40 


5 


37 


3.50 


2 


16 


2.015(5) 


6.74 


6 


35 


3.75 


2 


16 


2.014(6) 


1.68 


6 


95 


4.00 


2 


16 


2.011(7) 


2.75 


6 


84 


3.00 


3 


256 


3.011(16) 


0.66 


2 


72 


3.25 


3 


64 


3.012(11) 


3.18 


4 


53 


3.50 


3 


64 


3.014(13) 


3.67 


4 


45 


3.75 


3 


16 


3.034(9) 


2.32 


6 


89 


4.00 


3 


64 


2.988(19) 


1.66 


4 


80 


4.00 


4 


64 


4.002(25) 


3.01 


4 


56 



Table 21: Estimates for the limiting ratios Tint £/(l)/Tint,f'(A;) as L — >■ oo from fits to a 
constant A. Error bars are one standard deviation. 
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Fits to Tint.f = AL'' 


n h T J ■ 






1.25 1 [256] 


0.053(1) 1.27(1) 


3.76 1 5 


1.50 1 [256] 


0.104(2) 2.20(2) 


3.84 1 5 


I. lb i zoD 


O.ioo(zj 3.2o(^Oj 


U.zi i DO 


2.00 1 256 
2.00 2 [256] 


0.215(3) 4.47(9) 
0.224(2) 2.06(3) 


0.06 1 81 
6.91 1 0.9 


2.25 1 256 
2.25 2 [256] 


0.286(4) 5.35(14) 
0.287(3) 2.62(5) 


2.16 1 14 
10.68 1 0.1 


2.50 1 256 

Z.50 z z56 


0.353(6) 6.35(23) 

0.353(4) 3.15(8) 


0.02 1 90 

n n/^ 1 on 

0.06 1 80 


2.75 1 256 
2.75 2 [256] 


0.424(8) 7.19(34) 
0.421(5) 3.66(12) 


0.07 1 80 
3.02 1 8 


3.00 1 128 
3.00 2 128 

3.00 3 256 


0.505(6) 7.52(23) 
0.516(4) 3.55(8) 

0.505(6) 2.51(9) 


1.18 2 56 
1.46 2 48 
0.04 1 84 


3.25 1 128 

3.25 2 Z56 

3.25 3 256 


0.590(7) 7.66(30) 
0.5o5(^i0j 3.96(^23j 
0.585(8) 2.64(13) 


3.14 2 21 

0.003 1 96 

0.02 1 90 


3.50 1 128 

3.50 2 16 
3.50 3 64 


0.676(10) 7.77(39) 

0.693(2) 3.54(2) 
0.695(3) 2.33(4) 


1.07 2 58 
2.39 5 79 
0.05 3 100 


3.75 1 32 
3.75 2 32 
3.75 3 128 


0.779(4) 7.14(13) 
0.784(3) 3.46(4) 
0.798(7) 2.14(8) 


0.94 4 92 
1.63 4 80 
0.70 2 71 


4.00 1 128 

4.00 2 128 
4.00 3 64 
4.00 4 256 


0.916(16) 5.44(46) 
0.909(12) 2.84(17) 
0.890(5) 2.09(5) 
0.938(17) 1.18(12) 


0.15 2 93 

0.75 2 69 
1.58 3 66 
0.06 1 81 



Table 22: Estimates for 2;int,f' from fits to Tjnt.f = AL^ . Error bars are one standard 
deviation. Lmin in brackets indicates a poor fit (confidence level < 10%). 
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Fits to Ti„t,£' = AL^ + B 


n h T J ■ 


z- f> A B 


DF CL(%) 


1.25 1 128 


-0.213(35) -1.61(5) 2.19(8) 


0.49 1 48 


1.50 1 32 


-0.034(7) -15.86(2.92) 17.09(2.96) 


0.82 3 85 


1 TtX 1 1 OO 

i./o i izs 


O.Oool^zoj 9.D6(4.4zj — /.o6(4.ozj 


n OQ 1 c:t 


2.00 1 16 
2.00 2 32 


0.145(4) 10.37(51) -8.44(57) 
0.141(5) 5.48(32) -4.83(36) 


2.32 4 68 
1.32 3 72 


Z.Zo 1 32 
2.25 2 [128] 


0.235(8) 8.92(61) —6.68(86) 
0.200(19) 6.37(1.31) -6.45(1.89) 


1.62 3 65 
3.44 1 6 


2.50 1 32 
2.50 z 64 


0.315(8) 9.01(59) -6.75(1.00) 

0.294(11) 5.38(48) —5.23(85) 


0.21 3 98 

0.91 2 63 


2.75 1 16 
2.75 2 128 


0.411(5) 8.19(30) -4.76(57) 
0.348(24) 6.87(1.40) -9.62(3.28) 


1.06 4 90 
0.47 1 49 


3.00 1 32 
3.00 2 16 
3.00 3 16 


0.481(10) 9.08(63) —6.79(1.68) 
0.500(4) 4.00(11) -2.15(25) 

0.497(3) 2.72(6) -1.71(13) 


0.84 3 84 
1.35 4 85 
2.40 4 66 


3.25 1 64 

O.Zo Z oZ 

3.25 3 64 


0.558(21) 9.75(1.44) -11.92(5.96) 

n txw/'o\ A 0/i/'oo\ o ciof^o\ 
0.5M(oj 4.24(^22j — 2.93(/oj 

0.565(12) 3.10(26) -3.51(1.12) 


2.93 2 23 
0. ( i o o ^ 
0.22 2 90 


3.50 1 64 
3.50 2 16 
3.50 3 32 


0.648(25) 9.50(1.56) -14.38(8.87) 
0.691(6) 3.57(11) -0.13(40) 
0.689(8) 2.42(11) -0.66(55) 


1.15 2 56 

2.27 4 69 
0.34 3 95 


3.75 1 16 
3.75 2 16 
3.75 3 16 


0.790(9) 6.69(32) 2.96(1.49) 
0.800(6) 3.15(10) 2.21(50) 
0.795(5) 2.15(6) 1.10(28) 


0.33 4 99 
1.62 4 81 
2.91 4 57 


4.00 1 32 
4.00 2 32 
4.00 3 16 
4.00 4 128 


0.935(18) 4.76(48) 20.52(5.09) 
0.925(13) 2.54(18) 7.49(1.85) 
0.902(6) 1.94(6) 2.60(36) 
0.993(36) 0.80(19) 18.49(6.33) 


0.66 3 88 
0.59 3 90 
2.07 4 72 
0.07 1 79 



Table 23: Estimates for Zint,£' from fits to Tint,^' = AL^ + B. Error bars are one standard 
deviation. Lmin in brackets indicates a poor fit (confidence level < 10%). 
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Fits to Tint,f' = AL"" 


Fits to 






a/u 


f3/u 


Q 


-^min 


z 




z 


Difference 


(exact) 


(exact) 


1.25 






128 


-0.213(35) 




-0.35527 


0.11118 


1.50 






32 


-0.034(7) 




-0.22663 


0.11678 


1.75 


256 


0.158(2) 


128 


0.085(26) 


0.074 


-0.10929 


0.12131 


2.00 


256 


0.215(3) 


16-32 


0.143(3) 


0.072 


0.00000 


0.12500 


2.25 


256 


0.286(4) 


32 


0.235(8) 


0.051 


0.10363 


0.12798 


2.50 


256 


0.353(3) 


32-64 


0.307(7) 


0.046 


0.20357 


0.13034 


2.75 


256 


0.424(8) 


16-128 


0.408(5) 


0.016 


0.30168 


0.13212 


3.00 


128-256 


0.511(3) 


16-32 


0.497(3) 


0.013 


0.40000 


0.13333 


3.25 


128-256 


0.587(5) 


32-64 


0.572(7) 


0.015 


0.50126 


0.13393 


3.50 


16-128 


0.693(1) 


16-64 


0.689(4) 


0.004 


0.61007 


0.13377 


3.75 


32-128 


0.784(2) 


16 


0.796(4) 


-0.012 


0.73760 


0.13242 


4.00 


64-256 


0.898(4) 


16-128 


0.910(5) 


-0.012 


1.00000 


0.12500 



Table 24: Summary of estimates for -2int,f', averaged over all values of k for which the 
confidence level is > 10%. Error bars are one standard deviation. 
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t/Ch 


= ^ log L + 


B 


T 


ICh=ALp 




q 


k 


-^min 


A 


B 


CL(%) 


P 


A 


CL(%) 


1.75 


1 


256 


0.159(8) 


2.28(5) 


9 


0.049(3) 


2.41(4) 


14 


2.00 


1 


128 


0.224(8) 


2.74(4) 


23 


0.056(2) 


2.91(3) 


44 


2.00 


2 


128 


0.125(3) 


1.23(1) 


8 


0.064(1) 


1.34(1) 


30 


2.00 


2 


256 


0.134(5) 


1.18(3) 


59 


0.067(2) 


1.32(2) 


36 


2.25 


1 


128 


0.337(13) 


3.11(7) 


22 


0.068(3) 


3.42(5) 


45 


2.25 


1 


256 


0.373(25) 


2.89(14) 


86 


0.072(5) 


3.33(9) 


69 


2.25 


2 


128 


0.178(5) 


1.45(2) 


30 


0.073(2) 


1.63(2) 


25 


2.50 


1 


64 


0.410(13) 


3.86(6) 


33 


0.069(2) 


4.19(4) 


79 


2.50 


1 


128 


0.438(22) 


3.71(12) 


62 


0.071(3) 


4.13(8) 


85 


2.50 


2 


64 


0.215(5) 


1.85(2) 


15 


0.073(2) 


2.02(2) 


76 


2.50 


2 


256 


0.244(14) 


1.68(8) 


55 


0.077(4) 


1.97(5) 


74 


2.75 


1 


64 


0.531(21) 


4.55(10) 


30 


0.073(3) 


4.99(7) 


65 


2.75 


1 


128 


0.584(35) 


4.26(19) 


89 


0.078(5) 


4.87(12) 


88 


2.75 


2 


64 


0.286(7) 


2.15(3) 


62 


0.079(2) 


2.41(2) 


36 


3.00 


1 


32 


0.619(20) 


5.63(9) 


61 


0.073(2) 


6.05(6) 


89 


3.00 


1 


64 


0.660(32) 


5.42(15) 


100 


0.075(4) 


5.98(10) 


90 


3.00 


2 


64 


0.354(11) 


2.57(5) 


32 


0.081(2) 


2.89(4) 


82 


3.00 


2 


128 


0.383(19) 


2.41(10) 


94 


0.084(4) 


2.83(6) 


98 


3.00 


3 


64 


0.240(6) 


1.68(3) 


28 


0.083(2) 


1.90(2) 


80 


3.00 


3 


128 


0.257(11) 


1.59(6) 


98 


0.085(3) 


1.88(3) 


82 


3.25 


1 


32 


0.833(30) 


6.30(13) 


22 


0.083(3) 


6.91(9) 


41 


3.25 


1 


64 


0.907(48) 


5.93(23) 


55 


0.086(4) 


6.79(15) 


38 


3.25 


2 


32 


0.421(11) 


3.12(4) 


48 


0.084(2) 


3.43(3) 


99 


3.25 


2 


64 


0.443(17) 


3.01(8) 


81 


0.084(3) 


3.42(5) 


98 


3.25 


3 


64 


0.311(9) 


1.92(4) 


85 


0.089(3) 
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0.124(5) 


2.41(5) 


63 



Table 25: Fits for the ratio ri^t,e' /Cr- The quoted error bar corresponds to one standard 
deviation. 
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Fits to rint,C2 = 


n h T J ■ 




DF CL(%) 


1.25 1 256 


0.010(1) 1.05(1) 


0.05 1 82 


1.50 1 [256] 
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2.00 1 256 
2.00 2 [256] 
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0.897(16) 2.99(25) 
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0.96 2 62 
1.25 2 54 
0.04 1 85 



Table 26: Estimates for z^tfi^ from fits to Tint,C2 = AL^. Error bars are one standard 
deviation. Lmin in brackets indicates a poor fit (confidence level < 10%). 
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Table 27: Estimates for Zint,C2 from fits to rint,C2 = ^-^^ + Error bars are one standard 
deviation. Lmin in brackets indicates a poor fit (confidence level < 10%). 
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Table 28: Summary of estimates for Zi^t^^, averaged over all values of k for ■which the 
confidence level is > 10%. Error bars are one standard deviation. The last column of 
each group sho'ws 2;int,C2 ~ ^int,£' from the corresponding fits. 
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Figure 1: Fits of ^/L to x* + BL-p. (a) q = 1.25, (b) g = 2, (c) g = 2.75, (d) 
q = 3.25, (f) q = 3.5. For q = 3.25 the fit is actually to a constant x*. 
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Figure 2: Fits of x/L^'" to A + BL^^. (a) q = 1.25, (b) g = 2, (c) g = 2.25, (d) g = 3. 
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Figure 3: Normalized autocorrelation function pooif) for the observables O = S', T' and 
C2, in the case g = 3, A; = 3, L = 256. 
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Figure 4: Our best estimates of the dynamic critical exponent ^int,^'- The purely statistical 
error bar is indicated in black and is usually about the same size as the symbol; the 
combined statistical and systematic error bar is indicated in blue. (For q = 1.25,1.50 
we are unable to estimate the systematic errors.) The static exponents a/u and P/u are 
shown for comparison (red and green curves, respectively). 
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Figure 5: Fits to ri^t.^ = AL^/" + B + CL-p for g = 2 (A; = 1). (a) 
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